The Scott conjecture for large Coulomb systems: a review

We review some older and more recent results concerning the energy and particle distribution in ground states of heavy Coulomb systems. The reviewed results are asymptotic in nature: they describe properties of many-particle systems in the limit of a large number of particles. Particular emphasis is put on models that take relativistic kinematics into account. While non-relativistic models are typically rather well understood, this is generally not the case for relativistic ones and leads to a variety of open questions.


Introduction and historical background
1.1. Many-particle quantum mechanics. Properties of ground states of large Coulomb systems involving N electrons, such as atoms or molecules, are of fundamental interest in quantum physics and chemistry. Notable examples are the ground state energy and the electron distribution in the ground state. The latter may be expressed in terms of the one-particle ground state density, i.e., the probability density of finding one of the N electrons at a specific location in R 3 . It is well known that systems on atomic length scales are accurately described by quantum mechanics [105,106]. This understanding relies on precise investigations of the underlying Hamilton operator.
We consider a molecule that consists of K point-like nuclei of charges Z = (Z 1 , ..., Z K ) ∈ (0, ∞) K , fixed at pairwise different positions R = (R 1 , ..., R K ) ∈ R 3K , as well as N electrons, all interacting via Coulomb potentials in the Born-Oppenheimer approximation. The total nuclear charge is |Z| := K κ=1 Z κ . The number of spin degrees of freedom is denoted by q ∈ N. Although in reality q = 2, one may, for notational convenience, choose q = 1 when the spin-dependence is trivial.
A non-relativistic quantum mechanical description of this system is provided by the operator We choose Hartree units, so that = e = m = 1, where , e, and m, denote the rationalized Planck constant, the elementary charge, and the electron mass, respectively. In the atomic case (K = 1, R = 0), we have U = 0 and write H N,Z ≡ H N,Z/|x| . Since electrons are fermions, they obey the Pauli exclusion principle, i.e., the Hilbert space in which the operator (1.1) acts is given by N ν=1 L 2 (R 3 : C q ), i.e., the subspace of L 2 (R 3N : C q N ) consisting of all square-integrable, C q N -valued functions whose sign changes under the exchange of any two particle coordinates.
We write for the lowest spectral point of the Hamiltonian H N,V . This number E S (N, Z, R) may or may not be an eigenvalue, and, if it is, it may be degenerate. While the results in this review concern (1.4), there is an important, related quantity, which has not received the mathematical attention it deserves; see (5.1) below.
In addition to ground state energies, we will be interested in one-particle ground state densities. We recall that the one-particle density of a general (pure) state for its density and analogously for other Hamiltonians that we discuss later. Although this is might be an abuse of notation since the eigenvalue E S (N, Z, R) could be degenerate, our statements about ρ S will be true for any choice of an eigenfunction. The notion of a one-particle density and, in particular, of a one-particle ground state density can be generalized to the case of mixed states, but we do not do this here. Also, if the lowest point in the spectrum E S (N, Z, R) is not an eigenvalue, one can still obtain meaningful statements for so-called approximate ground states, but we will not discuss this in this introduction. The goal of this review is to summarize known results and open questions concerning the ground state energy and the one-particle ground state density in the limit of large electron numbers and nuclear charges for non-relativistic and, especially, for certain relativistic descriptions of Coulomb systems. In the rest of this introduction we will focus on results for non-relativistic atoms. This and other settings will be treated in more detail in later sections, see the table of contents and Subsection 1.6 for relevant pointers. (1) It is well known that the spectral analysis of N -particle systems for fixed N is prohibitively difficult already when N ≥ 2, since the O(N 2 ) many interparticle interactions prohibit a reduction to a three-dimensional (possibly) soluble one-particle problem. (For instance, if the electron-electron repulsion was absent and K = 1 in H N,V , then one could separate variables to end up with the direct sum of Schrödinger operators describing hydrogen.) Instead, one often considers the properties of a system for a large number of particles. This leads to the study of asymptotic properties. In this review we entirely focus on results in the limit Z 1 , ..., Z K , N → ∞. The precise way of carrying out this limit when K > 1 is explained later. (2) Studying asymptotics clearly leads to less quantitative mathematical statements and is also questionable from a physical point of view since experimentally observed values of Z are bounded, e.g., by 92 for stable atoms. However, the mathematical analysis is drastically simpler and, interestingly, leads to theorems that coincide astonishingly well with experimentally measured data. (This observation has been made repeatedly in different contexts in mathematical physics. Stell [227, p. 48] calls it the principle of unreasonable utility of asymptotic expansions and makes some interesting philosophical remarks.) (3) There are some notable exceptions, however.
(a) For instance, recently much progress has been made in the investigation of smoothness properties of single eigenfunctions and sums of squares of eigenfunctions of many-particle Coulomb Hamiltonians, such as H N,Z , for fixed N . In this regard see, e.g., the works [75,78,76,77,81] for non-relativistic and [80] for (pseudorelativistic) Chandrasekhar atoms. Such a priori estimates for many-particle eigenfunctions are important, e.g., for the derivation of eigenvalue asymptotics for the associated one-particle density matrix [220]. (b) Another example concerns the maximal ionization of atoms (and molecules). Experiments indicate that doubly or higher charged anions do not exist (Massey [164,165]), i.e., one expects at most Z + 1 many electrons to be bound to the nucleus, while any further electrons are located infinitely far away with vanishingly small kinetic energy. Proving this claim is a notoriously difficult problem in mathematical physics, see, e.g., Nam [179,180] for recent reviews. A slightly weaker formulation, the so-called ionization conjecture, states that there is a number Q < ∞ such that, if E S (N, Z) is an eigenvalue, then N ≤ Z + Q. Two well known results in this direction are due to Lieb [145] and Fefferman and Seco [66,67], who proved that N < 2Z +1 and N ≤ Z + CZ 47 56 , respectively, are necessary conditions for E S (N, Z) to be an eigenvalue. Recently, Nam [178] improved Lieb's result and showed N < 1.22Z + 3Z 1 3 , which leads to a sharper result when Z ≥ 6. Lieb's result implies the fact that doubly negatively charged hydrogen atoms do not exist.
1.2. Glimpse at Thomas-Fermi density functional theories. The N particle quantum Coulomb problem of computing E S (N, Z, R) and the associated eigenspace is -like its classical analogue, the Kepler problem -prohibitively difficult to solve (even numerically) already for N ≥ 2 because of the O(N 2 ) many interactions between the N electrons. This necessitates the derivation of so-called "effective theories", i.e., energy functionals or equations, which depend only on a fixed, but small number of variables, like three or six, and describe at least the macroscopically observed properties of the given system "sufficiently accurately". Although these theories are usually more accessible to numerical analysis, they also pose some interesting mathematical challenges in view of the presence of nonlinearities, which simulate the interparticle interactions. Here we focus on so-called density functional theories, i.e., energy functionals, that only depend on the oneparticle density.
For simplicity, assume from now on the neutral, atomic case (K = 1, N = Z). The breakthrough in the description of ground state properties of H N,V came with the help of a particularly simple density functional theory, the so-called Thomas-Fermi theory [234,71,72], which will be reviewed in Subsection 2.1. In their seminal work [152], Lieb and Simon connected Thomas-Fermi theory to the quantum problem of finding E S (Z) and showed that the Thomas-Fermi energy E TF (Z), i.e., the infimum of the Thomas-Fermi functional E TF Z (Lenz [138]), is the leading term of the asymptotic expansion of E S (Z) when Z → ∞. The Thomas-Fermi energy scales like E TF (Z) = E TF (1) · Z 7/3 , which is a consequence of . Thus, the result of Lieb and Simon for the ground state energy reads E S (Z) = E TF (1) · Z 7/3 + o(Z 7/3 ) as Z → ∞, (1.7) see also Theorem 3.1. A numerical computation shows that E TF (1) ≈ −0.484 29 · q 2/3 , cf. Gombás [94, p. 60]. Figuratively speaking, the leading order in (1.7) is generated by the bulk of the electrons, which are located on distances O(Z −1/3 ) from the nucleus, and are described semiclassically. It should not come as a surprise that this energetic result is accompanied by a result connecting the quantum ground state density ρ S with the minimizer of E TF Z , the Thomas-Fermi density ρ TF Z . Indeed, Lieb and Simon [152], and Baumgartner [12] showed that the suitably rescaled ground state density ρ S converges to the minimizer of the Thomas-Fermi theory for hydrogen. More precisely, one has, due to the scaling properties of Thomas-Fermi theory, the convergence lim Z→∞ Z −2 ρ S (Z −1/3 · ) = ρ TF Z=1 (1.8) when both sides are integrated against characteristic functions of bounded, measurable subsets of R 3 . In the context of the ionization conjecture, Fefferman and Seco [66] obtained (as a corollary) the convergence in a stronger topology, namely in the so-called Coulomb norm; see (2.6). The precise result is contained in Theorem 3.1.

1.3.
Quantum effects close to the nucleus. Although Thomas-Fermi theory correctly predicts the leading order of E S (Z), it turned out that, as Scott [200, p. 859] wrote in 1952, the Thomas-Fermi energy gives values for the binding energy (−1) · E S (Z) "which are too high by roughly 20%. The actual binding energies increase quite smoothly with increasing Z, which suggests the existence of a more appropriate formula." Naturally, this defect of simple Thomas-Fermi theory triggered some discussions. One year before Scott's publication, Foldy [74] had proposed the formula E S (Z) = c 1 · Z 12/5 + c 2 (Z). Here c 2 (Z) depends on E S (2) and the sum of the ionization potentials of all atoms with atomic number greater or equal three and less or equal Z. (Foldy does not give a bound on c 2 (Z) but seems to assume that c 2 (Z) = o(Z 12/5 ).) More importantly for us, c 1 is a constant that only depends on the chosen units and obeys c 1 > E TF (1). The exponent 12/5 was derived from numerical values of the electrostatic potential close to the nucleus as a function of Z (Dickinson [35]). Since these values were only available for Z ≤ 80, Foldy's formula was not expected to hold asymptotically as Z → ∞. In the discussion of his formula Foldy [74, p. 398] points out that Thomas-Fermi theory does not correctly take into account the following two effects close to the nucleus. On the one hand, such electrons are bound stronger to the nucleus, but, on the other hand, they screen the bulk of the electrons at larger distances to the nucleus. Foldy suspected the screening to dominate, which explains the inequality c 1 > E TF (1) despite the fact that Z 12/5 ≫ Z 7/3 for Z ≫ 1.
Scott [200, p. 867] made Foldy's observations more precise and suggested a different formula for E S (Z). He believed that Thomas-Fermi theory does correctly describe the leading order of the ground state energy expansion, but that the failure of Thomas-Fermi theory "is due partly to the shortcomings of the statistical model in the region nearest the nucleus, and partly to the effect of exchange". Like Foldy, Scott suggested that the few, but high-energy electrons that are located close to the nucleus should generate this correction. Due to their proximity to the nucleus, these electrons should be described quantum mechanically. Since the correction would be generated only by "finitely many" electrons, the electron-electron repulsion should be irrelevant and the order of the correction should be O(Z 2 ), i.e., in agreement with the magnitude of the eigenvalues of the hydrogenic Hamiltonian If one drops the electron-electron repulsion in the Hamiltonian (1.1), the corresponding ground state energy will also behave to leading order like a constant times Z 7/3 as Z → ∞, but with a constant different from E TF (1). There will also be a subleading correction, given by a constant times Z 2 , and, remarkably, the constant here is the same q/4 as in (1.10). This is not a coincidence and will become clear in the discussion below. In the same decade, Hughes [113,114] (lower bound) and the authors of [208,209,210] (upper and lower bound) proved this conjecture. That is, they rigorously derived the expansion (1.10); see Theorem 3.4. The proof in [208,209,210] relied in part on the mathematical and physical intuition gained in the precursor [213], where the Scott correction is proved in the absence of electron-electron repulsion. We will present this motivating result and its short proof in Subsection 2.3. Remark 1.3. As has been observed, e.g., by Conlon [28], Huxtable [115], or Sobolev [219], the Z 2 -correction is a consequence of the singularity of the Coulomb potential and cannot be explained semiclassically. For instance, Huxtable's result [115,Theorem 7] as Z → ∞ for any potential W ∈ C ∞ (R 3 ) satisfying c|x| 2 ≤ −W (x) ≤ C|x| 2 and |∇W (x)| ≤ c ′ |x| for some c, C, c ′ > 0. Here c TF is related to the non-relativistic Thomas-Fermi theory with potential W .
1.3.2. Strong Scott conjecture. The Scott conjecture has a close relative that concerns the ground state density, the so-called "strong form of the Scott correction" (in short "strong Scott conjecture" from now on). It was formulated by Lieb [144, pp. 623-624]; see also Heilmann and Lieb [104, p. 3629]. The strong Scott conjecture states that the suitably rescaled ground state density ρ S on distances of order Z −1 from the nucleus converges to q times the three-dimensional hydrogenic density, i.e., The latter is the sum of squares of the L 2 (R 3 : C)-normalized eigenfunctions ψ S n,ℓ,m of the hydrogen Hamiltonian The hydrogenic density ρ H S is rather well understood; see Theorem 3.11. In particular, the right side of (1.12) converges and is spherically symmetric. (The labeling of the eigenfunctions ψ S n,ℓ,m uses the decomposition into angular momentum channels and will be further explained in Subsection 3.1.2.) Note that S H is unitarily equivalent to Z −2 S H Z by scaling x → x/Z, where S H Z is defined in (1.9). Any eigenfunction ϕ Z of S H Z scales like ϕ Z (x) = Z 3/2 ϕ 1 (Zx), where ϕ 1 denotes the corresponding eigenfunction of S H 1 = S H . The 1996 work of Iantchenko et al. [117] showed, among other things, (1.14) see also Theorem 3.5. It follows from the convergence results there that the oneparticle ground state density ρ S is approximately spherically symmetric in the limit Z → ∞ on distances Z −1 from the nucleus.
1.4. Dirac-Schwinger correction. As mentioned in the previous subsection, Scott anticipated that the subleading terms in the expansion of the ground state energy E S (Z) should take into account the extreme quantum effects close to the nucleus, but also the exchange energy of the electrons, as proposed by Dirac [39]. On a formal level, Schwinger [199], as well as Englert and Schwinger [45,46,47] (see also Englert [44] for a textbook treatment) derived the third term in the asymptotic expansion of the ground state energy, which grows like Z 5/3 . In fact, this term is not only generated by the exchange energy of the electrons, but is also due to the semiclassical asymptotics of the eigenvalue sum of the operator − 1 2 ∆− Φ TF with the semiclassical parameter Z −1/3 and the Z-dependent Thomas-Fermi potential Φ TF (see (2.15)). In view of the results in Subsection 1.2, the occurrence of − 1 2 ∆ − Φ TF in the analysis of E S (Z) is not unexpected. A decade later, Schwinger's and Englert's derivation was made mathematically rigorous in the monumental work of Fefferman and Seco [68,61,70,64,62,63,65]; see also Bach [3,4] and Graf and Solovej [97] for simplifications and improvements of parts of Fefferman's and Seco's arguments and see [60] for a review of Fefferman's and Seco's proof. They proved the existence of a constant C DS > 0, which can be computed in terms of the Thomas-Fermi density ρ TF Z=1 , see [68, p. 528], such that In [60, pp. 6, 9-10] and [64, pp. 13-14] Fefferman and Seco make a conjecture concerning a fourth, possibly oscillating term in the expansion of E S (Z). Córdoba et al. [29,30,31] analyzed this term in detail and showed, in particular, that it is bounded from below and above by constants times Z 3/2 . 1.5. The necessity of a relativistic description. From a physical point of view it is questionable whether one can describe atoms with large nuclear charges nonrelativistically since already the bulk of the electrons is localized in orbitals whose distance to the nucleus is roughly Z −1/3 or less. As Z increases, the electrons become localized closer to the nucleus and, by Heisenberg's uncertainty principle, one expects that at least the velocities of the innermost electrons are a substantial fraction of the speed of light c. In fact, the non-relativistic energy for electrons on the length scale Z −1 in the field of a nucleus of charge Z is already −Z 2 /2. Thus, the virial theorem implies that its kinetic energy is Z 2 /2. In classical mechanics, this would show that the velocity of the electron is Z. Since the velocity of light is 137 in our units, a single electron in the field of a uranium nucleus (Z = 92) would therefore move with a speed of ≈ 92 137 · c, which indeed is a substantial fraction of the speed of light. For this reason, a relativistic description is mandatory, in particular on the short length scale Z −1 . Meanwhile, at distances Z −1/3 , electrons are expected to move with velocities 10% of the speed of light and, indeed, as we will see momentarily, relativistic effects are negligible to the leading, i.e., Thomas-Fermi order of E S (Z).
Before turning to details, we point out that already Scott [200, p. 866] anticipated possible shortcomings of his non-relativistic formula for large Z: "Relativity effects of all kinds have been disregarded so far. Though this simplification has no serious consequences for Z < 30, these effects are quite important for heavy elements. It would be a difficult task to calculate them accurately. A straightforward extension of Thomas' statistical method (Vallarta and Rosen [240]) is inapplicable to our present problem, because most of the correction originates in the region close to the nucleus where the statistical method is vitiated by the boundary effect, and, in fact, such methods would give an infinite binding energy. Moreover, the interaction between the electrons is not wholly electrostatic." Concerning an extension of "Thomas' statistical method " we refer to Subsection 2.4 for recent developments in this direction.
From a fundamental physical point of view, heavy atoms and molecules should be treated by relativistic quantum electrodynamics (QED) and the corresponding field theory. Unfortunately, many fundamental mathematical elements, e.g., the state space and the Hamiltonian, lack mathematical understanding. As a consequence, one is thrown back to approximate models. Here we will review three such approximate Hamiltonian models that have been derived by physical arguments from QED, and proven useful in applications. Moreover, we consider a mathematical simplification thereof and a density functional obtained in this vein.
The first model can be traced back at least to Chandrasekhar [24] in the context of stability of neutron stars. In it the single-particle kinetic energy −∆/2 is replaced by √ −c 2 ∆ + c 4 − c 2 with c being the velocity of light. Despite its mathematical simplicity, the resulting operator features many physical defects, such as the violation of the principle of locality. More crucially for us, it leads to ground state energies, that are much too low compared to experimental data and can only be applied to atoms with nuclear charge Z < 88.
Physically and chemically more accurate models are based on projected Coulomb-Dirac [37,38] operators, such as the Brown-Ravenhall [22] or the Furry [93] operator, which are applicable to atoms with nuclear charge Z < 125 and Z < 138, respectively. The latter is used in quantum chemistry to compute the ground state energy of large atoms or molecules to chemical accuracy, see, e.g., Reiher and Wolf [192].
A common property of relativistic operators is the fact that, at least for large momenta, the kinetic energy scales like the Coulomb potential. On a heuristic level, it is clear that the sole limit Z → ∞ is meaningless since the potential energy cannot be controlled by the kinetic energy anymore. Consequently, the total energy will not be bounded from below, and the atom becomes "unstable" for fractions Z/c beyond a critical model-dependent coupling constant. To make mathematically meaningful statements about asymptotics, one considers the limit when both Z and c tend to infinity simultaneously with a fixed ratio Z/c =: γ. (Of course, like the limit Z → ∞, the limit c → ∞ is questionable since c has a fixed value). The idea to introduce γ as a separate parameter goes back at least to Schwinger [198].
For γ ≤ 2/π, Sørensen [184] proved that in the above-described limit, the leading order of the ground state energy in the Chandrasekhar model is given by the Thomas-Fermi energy. Moreover, the ground state density on the Thomas-Fermi length scale converges weakly and in the so-called Coulomb norm (see (2.6)) to the hydrogenic Thomas-Fermi density [168,169]. This indicates that the bulk of the electrons on the length scale Z −1/3 does not behave relativistically.
On the other hand, electrons on the hydrogenic length scale Z −1 are located much closer to the nucleus and, as described before, are expected to lead to relativistic corrections of the Scott correction. In fact, Schwinger [198] derived a relativistic Z 2correction, which is lower than Scott's. This lowering was proved, by two different approaches, in Solovej et al. [223] and in [90]. Later, a relativistic correction of the Scott correction was also proved for the Brown-Ravenhall [91] and the Furry operator [100].
The relativistic generalization of the strong Scott conjecture was proved recently in [87] (see also [85]), i.e., the convergence of the suitably rescaled one-particle ground state density of Chandrasekhar atoms on the hydrogenic length scale 1/Z to the sum of the squares of the eigenfunctions of the one-particle Chandrasekhar operator. Shortly thereafter, the corresponding statement for the physically and chemically accurate Furry operator was proved [170]. These results underscores the fact that electrons close to the nucleus behave relativistically and that selfinteractions of the innermost electrons are negligible.
1.6. Organization. We briefly summarize the contents of the present review.
In Section 2 we review three examples of density functional theories. For the first two, we refer to March's and Lieb's reviews [161,144]; see also [206] for a recent review. First, and most important for us, we discuss Thomas-Fermi theory. Secondly, we review Weizsäcker's extension of Thomas-Fermi theory, which is physically and mathematically richer than Thomas-Fermi theory. Qualitatively, this extension correctly accounts for quantum effects of electrons close to the nucleus. Thirdly, we investigate the Hellmann-Weizsäcker functional, which served as the basis for Siedentop's and Weikard's proof of the Scott correction. Finally, we consider a density functional that reduces to the Thomas-Fermi-Weizsäcker functional in the non-relativistic limit. It was derived by Engel and Dreizler and was recently investigated from a mathematical point of view.
In Section 3 we consider non-relativistic atoms, ions, and molecules, both in the presence and absence of a self-generated field, and summarize theorems concerning the energy asymptotics and the convergence of the quantum density on both the Thomas-Fermi and the Scott length scales. Emphasis will be put on Scott's original derivation of the energy correction, as well as the initial motivating results in [213]. Section 4 is concerned with relativistic descriptions. We summarize results concerning the energetic asymptotics as well as the convergence of the density for all the three different relativistic models discussed in the introduction -the Chandrasekhar, the Brown-Ravenhall, and the Furry model.
In Section 5 we discuss some open questions.
Notation. We write A B for two non-negative quantities A, B ≥ 0 to indicate that there is a constant C > 0 such that A ≤ CB. If C = C τ depends on a parameter τ , we sometimes write A τ B. The notation A ∼ B means A B A. The indicator function of a set Ω is denoted by 1 Ω . The negative part of a real number or a self-adjoint operator A is defined by A − := max{0, −A} ≥ 0.

Density functional theories
In this section, we briefly review three examples of effective theories that are known to describe correctly at least the leading order of the ground state energy of large atoms and molecules. They are known as density functionals, i.e., energy functionals that only depend on the one-particle density of a given many-particle system. We refer to Lieb's detailed review [144] on Thomas-Fermi-type theories and to [206] for a recent survey of density (matrix) functional theories. In particular, [144] also treats extensions of Thomas-Fermi theory like Weizsäcker's inhomogeneity [239] and Dirac's exchange [39] correction. The first one will be of some interest for us since it generates a Scott correction, whereas the second one will be discussed in passing only.
2.1. Thomas-Fermi theory. We begin with the simplest non-relativistic "statistical model of the atom" (Fermi [71,72], Gombás [94]), which was formulated in the late 1920's independently by Thomas [234] and Fermi [71,72]. In the molecular case with K nuclei of charges Z = (Z 1 , ..., Z K ) ∈ (0, ∞) K situated at positions R = (R 1 , ..., R K ) ∈ R 3K , the so-called Thomas-Fermi (TF from now on) functional (Lenz [138]) is given by with V and U as in (1.2) and (1.3), respectively. In the atomic case K = 1, we write E TF Z ≡ E TF Z/|x| . The first term of E TF V (ρ) represents the kinetic energy and is derived via the following argument based on a semiclassical phase space integration. The TF model views the N non-relativistic quantum particles in a potential W as a classical gas in phase space. Since Planck's constant is h = 2π = 2π in our units, the density of the semiclassical gas is Thus, the semiclassical kinetic energy is with the Thomas-Fermi constant The second term in (2.1) represents the interaction energy between the electrons and the nuclei.
The third term in (2.1) is the electrostatic self-energy of the charge density ρ. It is defined, more generally, for ρ and σ by Note that (by Plancherel and the convolution theorem), D(ρ, σ) is sesquilinear and positive, and (by Cauchy-Schwarz), D(ρ, σ) ≤ D(ρ, ρ) · D(σ, σ). In fact, D(·, ·) defines a scalar product on the set I defined in (2.7) below. Thus, the right side of defines a norm on that space. This norm is sometimes called Coulomb norm. The TF functional (2.1) is defined on its natural domain (Simon [216]) i.e., for nonnegative densities with finite kinetic energy and finite electron-electron repulsion. These conditions automatically guarantee the finiteness of the electronnucleus-interaction. (The local singularities at the nuclei are controlled by the kinetic energy, whereas the long-range part is controlled by the electron-electron repulsion.) To describe a system of N electrons, we restrict the TF functional to the set Here, mathematically speaking, N need not be an integer.
In their seminal work [152], Lieb and Simon were the first to analyze this functional with mathematical rigor. (See also [144, Section II] and [206,Subsection 4.1] for more detailed reviews.) The following theorem asserts the existence and uniqueness of minimizers of the TF functionals on I and I N .
Then the following statements hold.
(1) (Unconstrained problem): There exists a unique 0 ≤ ρ TF (Z, R, x) such that (2) (Constrained problem): If N ≤ |Z|, then there exists a unique, non-negative is not a minimum, i.e., there are no negatively charged ions in TF theory.
(3) (Unconstrained Thomas-Fermi equation) In the unconstrained problem, the minimizer ρ TF ∈ I obeys ρ TF = |Z| and Moreover, if ρ ∈ I satisfies (2.11), then it minimizes E TF V on I. If K = 1, then ρ TF is spherically symmetric and decreasing. (4) (Thomas-Fermi equation in constrained problem) In the constrained problem with 0 < N ≤ |Z|, the minimizer ρ TF ∈ I N satisfies for some (unique) µ = µ(N ) ≥ 0. Moreover, there is no solution ρ ∈ I N to (2.12) for any µ other than ρ TF . When N = Z, then µ = 0, and otherwise µ > 0, i.e., E TF (N, Z, R) is strictly decreasing in N . As N varies from 0 to |Z|, µ varies continuously from ∞ to 0. Moreover, µ is a convex, decreasing function of N . (5) (Scaling) For any a > 0, the scaling relations hold. (3) If K = 1 in the unconstrained problem, then we write ρ TF (Z, 0, x) ≡ ρ TF Z (x) and E TF (Z, 0) ≡ E TF (Z). Similarly, in the constrained problem we shall write ρ TF (N, Z, 0, x) ≡ ρ TF Z (N, x) and E TF (N, Z, 0) ≡ E TF (N, Z). (4) The scaling relations for K = 1 show, in particular, that TF theory has "natural" length and energy scales, the Thomas-Fermi length scale Z −1/3 and the Thomas-Fermi energy scale Z 7/3 , respectively. For Z = 1 the minimizer ρ TF (either in the unconstrained or constrained problem) is called hydrogenic Thomas-Fermi density. The numerical value of the associated infimum is E TF (Z = 1) ≈ −3.678 74 · γ −1 TF , cf. Gombás [94, p. 60]. The following theorem due to Lieb and Simon [152,Theorem IV.5] (see also [144,Theorem 2.8]) summarizes some important properties of the TF density. (1) Let κ ∈ {1, ..., K} be arbitrary. Then, as x → R κ , one has In the ionic case (N < |Z|, µ > 0), ρ TF is compactly supported and C 1 away from the R κ . Then Φ TF obeys the Thomas-Fermi differential equation (2.16) (7) (Sommerfeld) In the neutral case (µ = 0) the Sommerfeld solution solves the TF differential equation (2.16) for |x| > 0 and x = R κ and it is the only power law that does so. Moreover, In the atomic case (K = 1), the TF density ρ TF obeys as |x| → ∞. In absence of electron repulsion, the Thomas-Fermi energy can be computed easily. This is important for the heuristic derivation of the Scott correction (Subsection 3.1.4) and is done in the following remark.
Remark 2.5 (Thomas-Fermi energy for the Bohr atom). Let K = 1 and Z > 0, and consider TF theory for an atom in absence of electron repulsion, i.e., For N > 0, let It is elementary to see that there is a unique minimizer ρ TF Bohr (N, Z) and that this minimizer satisfies the Euler-Lagrange equation with some µ > 0. Integrating the 3/2-th power of this identity leads to the relation , and then to the formula for the energy In 1935 Weizsäcker [239] proposed a correction of Thomas-Fermi theory that penalizes rapid changes of the density, which are expected to occur close to the nucleus.
Remark 2.6. Some words on the history: Weizsäcker introduced this correction to explain the rise of the mass defect per nucleon in a nucleus from very heavy (say uranium) to semi-heavy nuclei (like iron). To that end he consulted Gamow's liquid drop model for nuclei and argued that, as a consequence of the uncertainty principle, the "surface of the nucleus" must be smeared out. For otherwise, an instantaneous drop of the density with infinite slope would lead to an infinite kinetic energy, which is unreasonable. This smearing of the surface could be accounted for by replacing the eigenfunctions that are used in the derivation of the Thomas-Fermi functional, namely plane waves, by waves with linearly varying amplitude. This gives rise to Weizsäcker's term ρ −1 (∇ρ) 2 .
We consider the Thomas-Fermi-Weizsäcker (TFW from now on) functional . It is naturally defined on the set where the gradient is understood in the sense of distributions. For fixed particle number N ∈ (0, ∞), the functional is defined on Weizsäcker introduced (2.23) with A = 1. However, it is convenient to have A > 0 as an adjustable parameter, as we shall see soon.
The mathematical analysis of E TFW V started with the works of Benguria [13] and Benguria et al. [14]. Besides its mathematical richness, it turned out TFW theory describes -at least qualitatively -the physics of real atoms more closely than Thomas-Fermi theory. For instance, the TFW minimizer is finite at the nuclei and decays exponentially at infinity. Moreover, binding is possible and anions can be stable in TFW theory. For a concise summary of TFW theory we encourage the reader to consult [13,14], as well as [144, Sect. VII]. Here we restrict ourselves to a summary of the energy expansion and the minimizing density as |Z| → ∞. Our presentation follows closely Lieb [144] and Lieb and Liberman [146]. We start with the following result on existence and uniqueness of minimizers of the TFW functional.
(1) (Unconstrained problem) There is N c ∈ (|Z|, 2|Z|) such that the TFW functional E TFW V has a unique minimizer ρ TFW (Z, R, x) on A with particle number ρ TFW (Z, R, x) dx = N c . This minimizer satisfies the TFW equation The infimum is denoted by (2.28) (2) (Constrained problem) If N ≤ N c , then the TFW functional has a unique minimizer ρ TFW (N, Z, R, x) on A N . This minimizer satisfies the TFW equation with W (x) as in (2.27), µ ≥ 0, and µ = 0 for N = N c . The infimum is denoted by Remark 2.8. Benguria and Lieb [15] proved the previously-mentioned ionization conjecture for TFW molecules and showed 0 < N c − |Z| ≤ 270.74 · ( A 2γTF ) 3/2 · K. As we shall see below, the value A = 0.1859 is in some sense natural. Together with the value of γ TF = (6π 2 /2) 2/3 /2 this leads to the bound N c − |Z| < 0.7335 · K.
Theorem 2.7 shows, in particular, that anions can be stable in TFW theory. The next theorem says that ρ TFW on the TF length scale is described by ρ TF ; see also Solovej [221] for results when only some of the nuclear charges tend to infinity.
, and N > 0 so that λ := N/|Z| is fixed. Define z and r by Z = |Z|z and R = |Z| −1/3 r, respectively. Then weakly in L 1 if λ ≤ |Z| and weakly in L 1 loc if λ > |Z|. Naturally, the question arises of how close the two infima E TFW (N, Z, R) and E TF (N, Z, R) are, i.e., one seeks an upper bound on the right side of (2.32) Already in the neutral, atomic case (K = 1, N = Z) one might be tempted to say that the difference is O(Z 5/3 ) by plugging in the TF density ρ TF Z and using the scaling relation ρ TF However, this is not correct, as can be seen heuristically as follows. By Theorem 2.3, one has ρ TF 1 (x) ∼ const |x| −3/2 as |x| → 0, which makes it plausible that |∇ ρ TF 1 | ∼ const |x| −7/4 , but this is not square-integrable and leads to an infinite Weizsäcker term. Instead, as the following theorem shows, the difference E TFW (Z, Z, 0) − E TF (Z, Z, 0) is given, to leading order as Z → ∞, by a constant times Z 2 . This is the Scott correction in TFW theory. As in the quantum problem, the Z 2 -term originates from effects on the hydrogenic length scale Z −1 rather than the TF length scale Z −1/3 , see [144, p. 635]. Lieb [144,Theorem 7.30] also shows that the correction is independent of the electron number, i.e., it also holds when comparing the TF and TFW energies for ions with fixed ratio N/Z.
Let us return to the general, multi-center case. To describe TFW theory on the length scale Z −1 more precisely, we consider the atomic TFW functional without electron repulsion. After a 'renormalization' (that is, formally subtracting the integral of 2 5 γ TF Z γ TF |x| 5/2 from the right side of (2.23)) one can show that the resulting functional has a unique minimizer and that this minimizer solves the Euler-Lagrange equation see [144,Theorem 7.29]. By scaling, one has ρ( Then we have the following results on the hydrogenic energy and length scales. (2) (Density) Define z and r by Z = |Z|z and R = |Z| −1/3 r, respectively. Then the solution ρ TFW of the problem (2.30) converges to that of (2.33) in the sense that for each κ ∈ {1, ..., K}, both pointwise and in L 1 loc .
As emphasized by Lieb and Liberman in [146, Section 2.C], the second term in (2.34) has the following properties, which allows one to think of it as a "core effect": • It is independent of λ = N/|Z|. • It is additive in the nuclei, that is, it is a sum of terms corresponding to each atom in the molecule. • The constant D TFW does not change if the electron-electron repulsion is removed.
By the last point, we mean that (2.34) remains true, with the same constant D TFW , if in the definitions of both E TFW (N, Z, R) and E TF (N, Z, R) the term D(ρ, ρ) is dropped. This is proved in [144].
The asymptotics (2.34) and the convergence in (2.35) suggest a discussion of the parameter A in (2.23). By Theorem 2.10, it suffices to discuss the atomic case K = 1. While Weizsäcker initially chose A to be one, other values have been suggested. For instance, Kirzhnits [132] suggested A ≈ 1/9 based on the gradient expansion of the Hohenberg-Kohn functional, assuming the Coulomb potential was replaced by a "weak perturbing potential" ( [146, p. 12]). However, due to the local singularity, the Coulomb potential cannot be regarded as such a weak perturbing potential.
More then 15 years before Theorem 2.10 was proved, Yonei and Tomishima [249] analyzed (2.33) with µ such that the solution ρ ∞ obeys ρ ∞ = Z. From a numerical analysis they concluded A ≈ 1/5 (especially when Z > 25) leads to good agreement with the energy obtained from summing up the first Z eigenvalues of the hydrogen operator (1.9) (Bohr atom, cf. Remark 2.5 and Subsection 3.1.4).
Possibly inspired by Yonei's and Tomishima's work, Lieb and Liberman [146, (2.32)] chose A such that the Z 2 -correction in the TFW model agrees with that of the quantum model, i.e., D TFW = q/4. This choice leads to A = 0.1859.
Another choice for A is motivated by comparing the densities ρ TFW Z and the oneparticle ground state density ρ S in (1.5) on the length scale Z −1 . As indicated in (1.14), the spherical average over ρ S tends to the hydrogenic density ρ H S (cf. (1.12)) on the length scale Z −1 pointwise as Z → ∞. Recall that all hydrogenic eigenfunctions of (1.13) are finite at the origin with only eigenfunctions with ℓ = 0 being non-zero. (For a detailed analysis of ρ H S , we refer to Theorem 3.11 by Heilmann and Lieb [104].) Thus, the limiting value of Z −3 ρ S (Z −1 · (0+)) as Z → ∞ is well defined and can be computed explicitly thanks to the explicit knowledge of hydrogen eigenfunctions. At the same time, the convergence in (2.35) is pointwise as well, so Z −3 ρ TFW Z (Z −1 ·) is also accessible and can be computed numerically. Thus, to have agreement of the quantum density ρ S and the TFW density on the scale Z −1 , one may choose Weizsäcker's parameter A so that one has the equality In conclusion, one may regard the proof of the Scott conjecture in TFW theory as a warm-up problem for its proof in the full quantum problem. (One may wonder whether Scott was aware of Weizsäcker's extension [239] at the time of writing his work [200].) 2.3. Hellmann-Weizsäcker functional. In this subsection we discuss the socalled Hellmann-Weizsäcker functional, which plays an important role in the proof of the Scott conjecture by Siedentop and Weikard, as we will discuss in Subsubsection 3.1.4 below.
We first introduce the Hellmann-Weizsäcker functional for fixed angular momentum ℓ ∈ N 0 . We work on R + with the measure dr. We set Here, the derivative of √ ̺ ′ is understood in the sense of distributions and we recall that the square integrability of this derivative implies that √ ̺ is continuous on R + and has a boundary value √ ̺(0). In particular, the last condition in the definition of G W is well defined. Let α ℓ := π q(2ℓ + 1) 2 and define, for Z > 0 and ̺ ∈ G W , The second term is finite by Hardy's inequality and the last term is, since, for any Thus, E HW ℓ,Z is well defined on G W . Next, we introduce the full Hellmann-Weizsäcker functional. It is defined on The Hellmann-Weizsäcker functional [108] is defined, for ̺ ∈ M W , by One can prove that E HW Z is well defined on ̺ ∈ M W . This functional is studied in detail in [207]; see also Hoops [112]. Finally, for N > 0, we set The following theorem shows how the terms (α ℓ /3)̺ 3 ℓ in the Hellmann-Weizsäcker functional are related to the term (3/5)γ TF ρ 5/3 in the TF functional.
These functionals are well defined on sets G and M that are defined in a similar manner as G W and M W . A straightforward computation [213,Theorem 1] shows This is one step in the proof of Theorem 2.11.
(2) Hoops [112,Theorem 4.5] showed that the electron-electron repulsion does not alter (2.40) significantly. Moreover, he computed the coefficient of the Z 2 -term in this case. If Z − αZ β ≤ N ≤ Z + Q c , where α > 0, 0 ≤ β < 2/3, and Q c ≥ 0 is the Z-independent number specified in [112,Theorem 3.3], then one has where G is the infimum of another explicit functional defined in [112, (4.15)-(4.16)] and obeys the numerical bounds 2 · 0.388 ≤ G ≤ 2 · 0.417, see [112, p. 58]. The coefficient G is about three times bigger than Scott's coefficient 1/4. As Hoops puts it [112, p. 58]: "To have such a big discrepancy suggests that the Hellmann-Weizsäcker functional does not treat the innermost electrons sufficiently accurate to get the same behavior as the quantum mechanical ground state. That means that Weizsäcker's gradient term is a major correction (it creates a Z 2 -order term) at places where we have strong varying potentials but it does not suffice to give the right coefficient." 2.4. Relativistic TFW functional by Engel and Dreizler. As discussed in the introduction, a relativistic description of large Coulomb systems is mandatory. This suggests to consider relativistic density functionals. A particularly simple one can be traced back at least to Vallarta and Rosen [240] and Jensen [124], who mimicked the steps (2.2)-(2.4) with the kinetic energy p 2 /2 replaced by c 2 p 2 + c 4 − c 2 to derive a relativistic Thomas-Fermi theory. For q = 2 and nuclear configuration This functional is unbounded from below since the relativistic kinetic energy cannot control the Coulomb singularity. This was already anticipated by Jensen, see also Gombás [94, §14], [95, Chapter III, Section 16] for a review of these facts. Gombás also suggested that Weizsäcker's (non-relativistic) inhomogeneity correction would prevent the unboundedness from below. His suggestion was first carried out by Tomishima [237], who showed, among other things, the finiteness of the energy and the electron density at the nucleus. While Gombás introduced the Weizsäcker term ad hoc, Engel and Dreizler [43] offered a (formal) derivation from quantum electrodynamics. The Engel-Dreizler derivation also yields an exchange term. In total, their functional reads The Weizsäcker term is

and an adjustable parameter
The analysis of E rTFWD c,V started with Chen [25] and was continued in the works [27,26]. In the ultrarelativistic limit, i.e., in absence of the arsinh function in the Weizsäcker term, it had been investigated earlier in [16]. The functional E rTFWD c,V is naturally defined on In absence of the exchange term X (ρ), Chen [25, p. 39] proved the existence of minimizers of E rTFWD As we shall discuss in Subsection 4.1, non-renormalized relativistic quantum models for Coulomb systems are not expected to be well defined for arbitrary large nuclear charges, since the Coulomb potential and the kinetic energy have the same scaling behavior, at least for high momenta. The renormalization in Engel's and Dreizler's derivation leads to the arsinh function in Weizsäcker's term which ensures the lower boundedness of E rTFWD c,V (ρ) for all nuclear charges Z. The necessity of renormalization was realized early, see, e.g., Heisenberg and Euler [107]. . Let ξ := (4π) −1 max{X(t)/t 3 : t > 0} and s 0 : R + → R + be the explicit function given in [26, (16)], which is strictly monotone increasing and satisfies s 0 (0) = 0 and lim κ→∞ s 0 (κ) = ∞. Then for all ρ ∈ P with ρ = N one has In absence of the Dirac term, an analogous result was proved in [27,Theorem 1]. In fact, their result also holds in the molecular case.
The existence of minimizers and bounds for the excess charge (in absence or presence of the Dirac term) were proved by Chen [25, p. 39] and in [27,Theorem 2] and [26,Theorem 2]. In passing, we mention that bounds on the excess charge are available in many non-relativistic models; see, e.g., Lieb [144], as well as Benguria and Lieb [15], Solovej [221], and the more recent work [88].
The following result concerns the energy asymptotics in the atomic case.
Theorem 2.14 ([205, Theorem 1]). Let K = 1 and Z/c > 0 be fixed. Then The core elements of the proof are the facts that the relativistic kinetic energy is dominated by the non-relativistic one, that relativistically described electrons far away from the nucleus behave non-relativistically, and that the TFW functional provides a Z 2 -correction to the Thomas-Fermi energy. One may argue that Theorem 2.14 is expected in view of the heuristics explained in Subsection 4.3.1. Still, it is quite surprising that a formally derived functional correctly yields a fundamental feature like the ground state energy of large atoms.

Quantum mechanics of non-relativistic Coulomb systems
In this section we review results for the energy and density of non-relativistic systems with one or several nuclei. These are described by the Hamilton operator 3.1. Atoms without magnetic fields.
3.1.1. Thomas-Fermi scale. One of the main results of the work [152] of Lieb and Simon is that TF theory correctly describes both the leading order of the quantum mechanical ground state energy of large atoms, and the one-particle electron distribution on the length scale Z −1/3 in the limit Z → ∞. This is summarized in the following theorem.
is an eigenvalue of H N,Z and ρ S is the one-particle density of any of its associated normalized eigenfunctions, then one has In the context of a proof of stability of matter, Thirring [233] found a substantially simpler proof of the lower bound in (3.1) that used coherent states involving Gaussians. Lieb slightly generalized these coherent states and found a shorter proof for the upper bound in (3.1), see [144,Theorem 5.1] (where also an adaption of Thirring's proof of the lower bound is reported). Formula (3.3) (for certain characteristic functions) is also due to Baumgartner [12], while the convergence (3.2) of ρ S in Coulomb norm was proved by Fefferman and Seco [66]. (Note that the Scott correction with error bound O(Z 47/24 ) (see [208,209,1,2]) actually implies the quantitative bound follow from the fact that finite sums of characteristic functions of bounded measurable sets are dense in L 5/2 , together with the uniform boundedness (1), which in turn follows from the kinetic Lieb-Thirring inequality. . Let E TF Bohr (N, Z) be the corresponding quantity in TF theory, defined in Remark 2.5, and note that, by scaling, . Indeed, this can either be proved by following the proof of Theorem 3.1 or, much more directly, by using the explicit formula for the eigenvalue E Bohr (N, Z) and the explicit formula for E TF Bohr (N, Z) in (2.22). Let us return to the situation of Theorem 3.1. Recall that the TF minimizer satisfies ρ TF z=1 (x) ∼ const |x| −3/2 as x → 0 (see Theorem 2.3). Thus, Theorem 3.1 implies that the ground state density develops a singularity on the scale Z −1/3 . This is not surprising, since the electrons in a Bohr atom have a density proportional to Z 3 , which tends to infinity relative to the TF magnitude Z 2 . Thus, TF theory is not expected to provide a complete description of large atoms, especially on scales Z −1 close to the nucleus. For this reason we now look closer at the electrons on the shorter, hydrogenic length scale Z −1 .

Scott scale.
As explained in the introduction, Scott's idea is that the leading correction to TF theory stems from the few innermost high-energy electrons on distances Z −1 from the nucleus, i.e., on the natural length scale of the hydrogen operator S H Z . Recall also that the eigenvalues of this operator have magnitude O(Z 2 ).
The following theorem states that Scott's conjecture is indeed true. Hughes [113,114] (lower bound) and the works [208,209] (upper and lower bound) proved this conjecture for neutrally charged atoms. Since the Scott correction should not depend on "electrons on the outermost shells", it was believed that the conjecture also holds for ions. Indeed, Bach [2,1] showed that also this intuition is correct, and proved Scott's conjecture for (positive and negative) ions.
We will give a heuristic argument in favor of Theorem 3.4 later in Subsection 3.1.4. At this point we would also like to mention two alternative and interesting proofs of Theorem 3.4 by Ivrii and Sigal [119] and by Solovej and Spitzer [224], which do not use the spherical symmetry of the atomic case (see Theorem 3.15). Both used them to prove the Scott conjecture for clamped nuclei whose distances are scaled by Z −1/3 . Next, we turn to the strong Scott conjecture. Scott's heuristics led Lieb [144, p. 623] to the belief that for large Z the suitably rescaled electron density ρ S on distances Z −1 from the nucleus converges to the hydrogenic density ρ H S (see (1.12)). Since the magnitude of the total electronic density of a Bohr atom is proportional to Z 3 , one suspects that the correct object to study is the scaled density Z −3 ρ S (Z −1 x). A step towards the proof of the strong Scott conjecture in the full N -particle setting was taken in [202], where the following upper bound was shown: see also Theorem 3.8 for a more general statement. Interestingly, the numerical value π/24 is quite close to q −1 · ρ H S (0) = π −1 ζ(3) (see Theorem 3.11). Motivated by this strong numerical evidence, Iantchenko et al. [117] eventually proved Lieb's strong Scott conjecture.
In fact, in [117] also an angular-momentum-resolved version of this conjecture is proved. To formulate this result, for for ℓ ∈ N 0 and r > 0. (In the following we use the letter ̺ : R + → R + to denote one-dimensional, radial densities with particle number ∞ 0 ̺(r)dr, i.e., we integrate with respect to dr and not r 2 dr. Three-dimensional densities are denoted by the letter ρ : R 3 → R + .) The densities ̺ ℓ in (3.7) are related to the total density ρ in (1.5) by for (3.7) and analogously for other Hamiltonians that we discuss later.
The following theorem states the validity of the strong Scott conjecture.
Theorem 3.5 ( [117]). Let ℓ ∈ N 0 , E S (N, Z) be an eigenvalue of H N,Z , and ρ S and ̺ ℓ,S be the total and angular momentum resolved one-particle densities associated to any of its eigenfunctions, respectively. Then the following statements hold.
(1) (Convergence of angular momentum density). For all r > 0 one has pointwise convergence (2) (Convergence of total density) Let W ∈ L ∞ (S 2 ) and r > 0. Then the total density, when spherically averaged, converges pointwise to the hydrogenic density, i.e., Moreover, for any locally bounded v ∈ L 1 (R 3 ), one has Remark 3.6. Theorem 3.5 remains valid for the total and angular momentum resolved one-particle densities ρ S and ̺ ℓ,S of elements in any sequence of normalized functions ψ N , Such a sequence (ψ N ) N ∈N is sometimes called approximate ground state of H N,Z on the Scott scale.
Later, we shall see that the hydrogenic density obeys as |x| → ∞. Since Theorem 2.3 asserts that the TF density has the exact same asymptotic behavior, but for small |x| (on the scale Z − 1 3 ), Theorem 3.5 shows that there is a smooth transition of ρ S between the length scales Z −1 and Z − 1 3 .
Remarks 3.7 (More results on the quantum density). We summarize some further results for the many-particle ground state density.
(1) Recall that Theorem 3.5 does not say anything yet about the quantum density at the origin r = 0. Therefore, it is of interest to find quantitative upper bounds on ρ S (0). We record the following theorem.
Theorem 3.8 ( [202,187,188]). Let (ψ N ) N ∈N be an approximate ground state of H N,Z on the Scott scale in the sense of Remark 3.6 and ρ S be the oneparticle density of the element ψ N of that sequence. Fix any c > 0 and assume that N/Z > 0 is fixed or that N > Z − cZ 1/2 . Then one has For N = Z this was shown in [202] before the strong Scott conjecture [117] was proved. The result for ions is due to Rakowsky [187] and [188].
At the core of the proof of Theorem 3.8 lies a linear response argument and the inequality due to Hoffmann-Ostenhof et al. [110]. Theorem 3.8 supports the belief that the strong Scott conjecture holds for ions as well.
(2) The work [118] proved the convergence of the one-particle ground state density matrix on Scott's scale. Consider an approximate ground state  be the kernel of the associated approximate one-particle ground state density matrix γ S . Observe that 0 ≤ γ S ≤ q is a trace class operator. Recall the L 2 (R 3 )-normalized eigenfunctions ψ S n,ℓ,m of the hydrogen operator (1.13). Define the orthogonal projection onto the negative (purely discrete) spectral subspace of this operator by where the series is pointwise convergent (see Theorem 3.11 by Heilmann and Lieb [104]). Then the following result holds.
Theorem 3.9 ([118, Theorem 1]). Let K be a trace-class operator and γ S be an approximate one-particle ground state density matrix of the Schrödinger operator H N,Z in (1.1). Then one has the weak convergence (3.24) (3) As already indicated, there is a smooth transition of ρ S between the scales Z −1 and Z −1/3 (in view of Theorems 2.3 and 3.5, and Formula (3.19)). This may have led to Lieb's belief [144] that Z −3/2−3/(2δ) ρ S (x/Z δ ) should also behave like (1/γ TF ) 3/2 |x| −3/2 when δ ∈ (1/3, 1). In [116] Iantchenko succeeded in proving this conjecture using some of Ivrii's and Sigal's methods [119] adapted to the case of (neutral) atoms.  1.1)) on the Scott scale in the sense of Remark 3.6 and let ρ S be the one-particle density associated to ψ N . 0)) 3 and |ν| = ν 1 + ν 2 + ν 3 . We also use the notation x := 1 + |x| 2 for x ∈ R d .) Then for any δ ∈ (1/3, 1), one has It may very well be true that Theorem 3.10 also holds under less severe assumptions on the test function U . Theorem 3.10 can also be generalized to the molecular case, see [116,Theorem 7] for details. Further extensions were recently outlined by Ivrii [122,120,121].
Some elements of the proof of Theorem 3.5 will be given in the relativistic case, which is in spirit the same but technically more elaborate (Subsection 4.3.3). Here we just mention that all known proofs for the convergence of the density (as of this writing) rely on energetic results. It would be interesting to see if this scheme could be reversed.
Before we come to the arguments in favor of Theorem 3.4, we collect some properties of the hydrogenic density ρ H S . where the R n,ℓ are explicitly known [104, p. 3630]. Thus, by Unsöld's theorem, Making heavy use of properties of special functions, they carried out the ℓ-and n-summations (in that order), and proved the following theorem. (1) The series on the right side of (3.29) converges pointwise.
(2) One has the asymptotic expansion, as r → ∞, The first few coefficients are (3) Let λ ∈ (0, 1] and ρ TF (λ, 1, x) be the TF minimizer of E TF z=1 on I λ . Then for any λ ∈ (0, 1], is monotone decreasing and achieves its maximum at r = 0 with Remark 3.12. The "shell structure", i.e., the oscillations of ρ H S (r) are barely visible. As Heilmann and Lieb put it [104, p. 3633]: "In fact, it is necessary to take two derivatives with respect to r in order to make the oscillatory terms as large as the nonoscillatory ones. In short, shell structure is not a prominent property of this universal atomic function." For graphical illustrations, see, e.g., [104, p. 3633].

Scott's derivation of the Z 2 -correction and ideas of Siedentop and Weikard.
We first present Scott's [200] heuristic derivation of the Z 2 -correction to TF theory. For simplicity, we restrict ourselves to the neutral case N = Z. Our exposition follows March [161, pp. 8-11]. Afterwards we present the ideas in [208] to prove the upper bound (3.5).
Scott's idea is that the Z 2 -correction stems from the innermost electrons, which live roughly a distance Z −1 away from the nucleus. Since there are only "few" electrons on this scale, one expects the electron-electron repulsion to be irrelevant there, i.e., the electrons should be described by the hydrogen operator S H Z . These heuristics suggest that the first correction of the Thomas-Fermi energy for a large atom should be the difference between the "Bohr energy" (i.e., the sum over all hydrogen eigenvalues) and the Thomas-Fermi Bohr energy. We computed the latter in (2.22). In the neutral case, we obtain On the other hand, in the Bohr atom, each shell (indexed by n ∈ N) of energy −Z 2 /(2n 2 ) has qn 2 states. Thus, the N = Z many electrons occupy k shells, where k is determined by and where 0 ≤ ε < 1 is the fraction of the (k + 1)-st shell that is filled. One finds Thus, the energy of the Bohr atom is The second term here is exactly the Scott correction.
We now present some ideas from [208] that enter into the proof of the upper bound in (3.5). As usual, an upper bound will be derived using a suitable trial state. As we will describe later in more detail, this trial state contains, in addition to hydrogen orbitals, also so-called Macke orbitals. It is in connection with these Macke orbitals that the Hellmann-Weizsäcker functional defined in Subsection 2.3 comes into play [212].
To motivate why this functional is relevant, we recall that the TFW energy has an asymptotic expansion whose first term coincides with the TF energy and whose second term behaves like Z 2 . If one could show that the TFW functional (with some choice of A) provides an upper bound for the ground state energy E S (N, Z), then there might be a chance that this functional can be used to prove the upper bound in (3.5). Apart from a factor that decreases like N −2 and is irrelevant in the limit N → ∞, this was indeed proved in one spatial dimension by March [162] and March and Young [163]. They showed that the sum of the first N eigenvalues where A N is the one-dimensional analog of (2.25); see also [206,Section 3] Ladányi [136] and [212] for precursors) the following upper bound was proved in the atomic case: with an (unimportant) remainder term Here n ℓ = (q(2ℓ + 1)) −1 ∞ 0 ̺ ℓ (r) dr and ε ℓ := n ℓ − [n ℓ ] ∈ [0, 1). Physically, n ℓ has the interpretation as the number of electrons in the sector of angular momentum ℓ corresponding to fixed magnetic quantum number m and spin s. (In [207,213,208], this number is denoted by N ℓ,m,s .) The total number of electrons in the angular momentum channel ℓ is therefore q(2ℓ + 1)n ℓ . Note that, if n ℓ is integer, then R(̺) ≤ 0 and the term can be dropped. For general n ℓ , this term has the effect of slightly modifying the prefactor of the term involving Thus, we have obtained an upper bound of the correct order, but not yet with the correct coefficient. In order to obtain the correct coefficient, we need one more modification of the above strategy.
After this preparation we can give a brief outline of the proof in [208] of the upper bound in (3.5). The basic idea is to describe electrons close to the nucleus using hydrogen eigenfunctions, and electrons far from the nucleus using the Hellmann-Weizsäcker functional. The distinction between "close" and "far" electrons is implemented by an angular momentum cut-off L = [Z δ ] for an appropriate δ ∈ (0, 1). This is suggested by the solution of Kepler's problem, where the perihelion of a planet grows like the square of its angular momentum. For the hydrogen atom, this is reflected by the explicitly computable expectation values of powers of |x| in their eigenfunctions, see, e.g., Bethe [17, (3.19)-(3.27)].
In the following outline, we disregard the electron-electron repulsion. Treating this term requires some effort, but it does not affect the Scott correction and, therefore, ignoring this term helps to clarify the basic steps in the proof.
(1) Use the variational principle with a trial state that is made up of hydrogen eigenfunctions for angular momenta 0 ≤ ℓ < [Z  (3.5). We mention that Macke orbitals and the Hellmann-Weizsäcker functional can also be used to give a proof of the lower bound of (1.10); see [210] and, for an exposition of the basic ideas, [211]. Both in the proof of the upper and lower bound, the use of Macke orbitals is reminiscent of the use of coherent states by Lieb in his proof of the asymptotic exactness of TF theory.

Thomas-Fermi scale.
As in the atomic case, TF theory predicts the leading term of E S (N, Z, R) and the density ρ S on the TF length scale correctly. More precisely, fix R = (R 1 , ..., R K ) ∈ R 3K , Z = (Z 1 , ..., Z K ) ∈ (0, ∞) K , |Z| = K κ=1 Z κ , and the TF electron number λ. It is not necessary to assume λ ≤ |Z|. For each N = 1, 2, . . ., define a N = N/λ. In H N,V in (1.1) replace each Z κ by a N · Z κ and each R κ by a −1/3 N R κ . This means that the nuclei come together at the rate a −1/3 N ∼ |Z| −1/3 as N → ∞ when N/|Z| is kept constant. At this point we want to emphasize that this scaling of R κ , which has become customary in the mathematical literature, is motivated by mathematical considerations rather than by physical reality. However, neither the energy nor the density is expected to be close to the molecular ground state energy or ground state density. In fact the two energies are expected to differ already to leading order, since the minimal positions of the nuclei do not scale in this way.
In case E S (N, a N Z, a −1/3 N R) is an eigenvalue, let ψ N be an associated eigenfunction with corresponding one-particle density ρ S . If not, ρ S shall denote the one-particle density associated to the element ψ N of an approximate ground state (ψ N ) N ∈N in the sense of Remark 3.2.
By the scaling relations (2.13) for TF theory, we have We now make the connection between the full quantum problem (1.4) with frozen nuclear positions associated to H N,V , and TF theory by letting the total nuclear charge |Z| and the electron number N tend to infinity in such a way that the ionization degree N/|Z| is kept constant. The following theorem due to Lieb and Simon [152] generalizes Theorem 3.1. (1) The quantity a −7/3 N E S (N, a N z, a −1/3 N R) has a limit as N → ∞ and this limit coincides with E TF (λ, z, R).
(2) The scaled one-particle density a −2 N ρ S (a −1/3 N x) associated to a (possibly approximate) ground state on the TF scale has a limit as N → ∞. If λ ≤ |z|, then the convergence is weakly in L 1 and the limit is ρ TF |z| (λ, x). If λ > |z|, then the convergence is weakly in L 1 loc and the limit is ρ TF |z| (|z|, x). As Lieb [141, p. 560] notices: "Note that if λ > |z|, then this result says that the surplus charge moves off to infinity and the result is a neutral molecule. This means that large atoms or molecules cannot have a negative ionization proportional to the total nuclear charge; at best they can have a negative ionization which is a vanishingly small fraction of the total charge." Remark 3.14. Note that, if the R κ are kept fixed and unscaled, one ends up with isolated atoms in the limit N → ∞, see Lieb [141, pp. 559-560] for the precise statements.

Scott scale.
In addition to being independent of the absence or presence of the electron-electron repulsion and of the ionization degree, one expects the Z 2correction to be the sum of the Scott corrections of the atoms constituting the molecule, as long as the minimal internuclear distance is not too close to the Scott scale. In particular, this is expected for ground states of molecules where the interatomic distance is expected to be of order one in the scaling parameter of the nuclear charges.
The Scott conjecture for molecules with frozen nuclear positions was first proved by Ivrii and Sigal [119] using a multiscale analysis and microlocal techniques. Later, Solovej and Spitzer [224] found a proof that is partly similar to the multi-scale analysis in [119] using an interesting new coherent states method. Around the same time, Balodis [8] proved the Scott correction uniformly in K. The following theorem is the content of [224, Theorem 1].

3.3.
Molecules with self-generated magnetic fields. In this subsection we consider molecules in the presence of a classical magnetic field. Quoting Erdős and Solovej [53, p. 229]: "External magnetic fields were taken into account in [153,154] (homogeneous) and [52] (inhomogeneous), but subject to certain regularity conditions. Self-generated magnetic fields, obtained from Maxwell's equation are not known to satisfy these conditions." For this reason we shall consider self-generated magnetic fields here. (See (2) in Remark 3.16 below for an explanation of the word "self-generated".) To introduce the setting precisely, let where all derivatives are understood in the sense of distributions and where ∇ ⊗ A denotes the 3 × 3 matrix of all derivatives ∂ i A j . We set |∇ ⊗ A| 2 = 3 i,j=1 |∂ i A j | 2 and, with c > 0 denoting the velocity of light, consider the magnetic field energy The identity here is a consequence of the Coulomb gauge.
In the non-relativistic approximation, the one-particle kinetic energy is the magnetic Schrödinger or Pauli operator in L 2 (R 3 : C 2 ) depending on whether the particles are (effectively) considered spinless, or have spin-1 2 . Here σ = (σ 1 , σ 2 , σ 3 ) is the vector of Pauli matrices. The total energy of a non-relativistic molecule with charges Z = (Z 1 , ..., Z K ) ∈ (0, ∞) K fixed at positions R = (R 1 , ..., R K ) ∈ R 3K in a classical magnetic vector potential A ∈ A is then described by the Hamilton operator (L 2 (R 3 : C 2 )), (3.45) where V is the electron-nucleus interaction as in (1.2). In case the kinetic energy is described by the Pauli operator, we will need to impose an additional restriction on the quotient |Z|/c 2 . Indeed, if |Z|/c 2 is sufficiently small, then the quadratic form associated to (3.45) is bounded from below uniformly in A, see, e.g., [92,69,148,147]. Crucially, however, stability fails, if |Z|/c 2 is too large, see [157,48].
For both choices of the magnetic kinetic energy T (A) and each fixed A ∈ A, the operator H N,V,A is defined as the Friedrichs extension of the corresponding quadratic form defined on S(R 3 : C 2 ).
For fixed magnetic potential A ∈ A and fixed nuclear positions, the electronic ground state energy is The total ground state density with fixed nuclear positions arises from minimizing this energy with respect to A ∈ A, i.e., Erdős and Solovej [53] proved the following theorem.   (1) The theorem is independent of the existence and uniqueness (modulo gauge freedom) of the minimizer A. In fact, it is not clear whether the infimum is attained at A = 0 or at a non-trivial magnetic field.
(2) The threshold κ 0 for which the assertion of Theorem 3.18 is shown to hold is less than the number κ cr above which H N,V,A fails to be bounded from below.
(3) The theorem does not assert that S(κ) is strictly decreasing, although this is believed to be the case. In fact, it is conceivable that S(κ) is constant equal to 1/4 for all κ up to the critical value κ cr beyond which it is minus infinity.
Remark 3.20. Theorem 3.17 concerning the leading order holds also in the case where the A-field is quantized [53]. So far, the Scott correction for quantized Afields -both in the non-relativistic setting of Theorem 3.18 and in the relativistic one of Theorem 4.21 below -can only be proved with a low ultraviolet cutoff of the magnetic field which corresponds to a length scale that is longer than the Scott scale, i.e., would be of limited physical meaning.
We emphasize that all results in this subsection concern the energy. We are not aware of results concerning the density.

Relativistic Coulomb systems
In this section, we discuss relativistic models of large Coulomb systems and begin with an overview of the relevant underlying one-particle operators. 4.1. One-particle operators. We first introduce the relativistic one-particle operators that will later be used to construct the many-particle operators that we are mostly interested in. We state conditions on the coupling constants of the Coulomb potential for which the operators can be defined and recall some of their spectral properties. More detailed treatments are contained, e.g., in the textbooks by Balinsky and Evans [7] and Thaller [231], as well as in the paper by Matte and Stockmeyer [166] and the references therein.

Chandrasekhar operator.
Definition. The Chandrasekhar operator is the simplest relativistic operator discussed here. In the literature this operator is sometimes referred to as Herbst operator or pseudorelativistic operator. Its origins can be traced back at least to Chandrasekhar in the context of (in)stability of neutron stars [24] (see also [155,156]). The mathematical investigation of this operator started with the work of Herbst [109]; see also Weder [243] for electric potentials V (x) = |x| −β with β ∈ (0, 1). The operator is defined as the Friedrichs extension of the quadratic form -whenever it is bounded from below (see (4.3)) -associated to with form domain being the Schwartz space S(R 3 : C). In the following, we abbreviate the fractional Laplace by |p| := √ −∆.
Scaling. The operator has a natural length scale, namely c −1 . Indeed, scaling x → x/c and writing γ := Z/c shows that C H c,Z is unitarily equivalent to c 2 C H 1,γ =: Kato's inequality. The sharp Hardy-Kato-Herbst inequality -for short Kato's inequality -states that In fact, Raynal et al. [189] showed that the form is strictly greater than −1, even if γ = 2/π. Numerical evidence for this fact had been provided by Hardekopf and Sucher [101].
Domain considerations. The quadratic form domain of C H γ is H 1/2 (R 3 ) when γ < 2/π. For γ = 2/π the form domain is the closure of S(R 3 ) with respect to the norm ( u, C H γ u + u 2 ) 1/2 . In analogy to the local case, we believe that there are functions in the form domain of C H 2/π for which both sides of Kato's inequality are infinite and, therefore, that this form domain strictly contains H 1/2 (R 3 ). For domain considerations, see also Le Yaouanc et al. [137].
(2) Instead of comparing Bessel kernels of C H ℓ,γ and p ℓ , one can compare the Bessel kernels of C H ℓ,γ and p ℓ − γ/r. For ℓ = 0, the latter can immediately be derived using the spectral theorem and recent heat kernel bounds for |p| − γ/|x| by Bogdan et al. [18]. Since for ℓ ≥ 1 the Hardy potential is again an operator perturbation for γ < 1/2 (by Hardy's inequality), Bessel kernel bounds for p ℓ − γ/r are similar to those for p ℓ . We expect that this strategy allows to remove the ε > 0 in (4.18) when γ > (1 + √ 2)/4. It is an open problem to prove that the behavior of ̺ H ℓ,C and ρ H C at r = 0 for γ > (1 + √ 2)/4 is optimal.
(3) We do not expect Bessel kernel bounds for p ℓ − γ/r to yield precise (probably γ-dependent) bounds for ̺ H ℓ,C at the origin in the cases ℓ = 0 and γ ≤ 1/2, and ℓ ≥ 1 and γ ≤ 2/π, because the Bessel kernel bounds for p ℓ − γ/r are similar to those for p ℓ .
Ground state transform. An alternative representation of the operator (−∆) 1/2 − γ/|x| in L 2 (R 3 ) proceeds via the ground state transform. To that end, recall (4.15). The ground state transform makes use of the fact that the radial function x → |x| −σ is a (generalized) ground state for |p| − Φ(σ)|x| −α . It states that Note that for α = 2 the ground state transform has been known long before, see, e.g., [190, p. 169] for a textbook treatment when d = 3. For a further study of the ground state transform in the fractional case, see [89].
Spectrum. Although the eigenvalues λ Z,n,ℓ (n ∈ N 0 ) of C H c,Z are not explicitly known, the inequality p 2 + 1−1 ≤ p 2 /2 and the lower bound of [91,Theorem 2.2] imply the inequalities − Z 2 2(n + ℓ + 1) 2 ≥ λ Z,n,ℓ ≥ −const · Z 2 (n + ℓ + 1) 2 , (4.22) where the constant in the second inequality can be chosen independently of γ ∈ [0, 2/π]. The expression on the left side is just the n-th eigenvalue of the hydrogen operator (1.9) in angular momentum channel ℓ. Thus, although the relativistic eigenvalues are smaller than the non-relativistic ones, their magnitudes in Z are the same and many of their summability properties with respect to n are similar. Finally, the spectrum of C H γ in [0, ∞) is purely absolutely continuous, the singular continuous spectrum is empty, and there are no embedded eigenvalues [109,Theorem 2.3]. In particular, there is no zero eigenvalue, a fact that we will use later.
Physical shortfalls. Although C H γ is mathematically well understood, it has a number of physical deficits. For instance, the restriction γ ≤ 2/π implies that only atoms with nuclear charge < 88 can be described. Moreover, the predicted ground state energies for heavy atoms are much too low. This can already be anticipated by comparing the ground state energies for hydrogen with c = 1 and coupling γ close to 2/π, which are ≈ −0.5 in the Chandrasekhar model [189, p. 106] and ≈ −0.06 (cf. (4.35)) in the Dirac model, respectively.
Although the model is unsuitable for the quantitative description of systems with strong attractive external Coulomb forces, Chandrasekhar [24] used it successfully for attractive two-particle Coulomb forces in his Nobel prize winning estimate on the mass necessary to collapse a star to a white dwarf. Despite its mathematical simplicity and its success in correctly describing some qualitative features of relativistic Coulomb systems, it is desirable to examine models that also lead to quantitatively correct predictions for ground state properties. Such models are, e.g., based on the Coulomb-Dirac operator, which we discuss next.

Coulomb-Dirac operator.
Free Dirac operator. In 1928 Dirac [37,38] derived a Lorentz invariant equation of motion for quantum mechanical particles with spin moving in an external electromagnetic field, the so-called Dirac equation. We refer to [17,231] for comprehensive treatments. For a free particle, the equation reads with the Dirac matrices α = (α 1 , α 2 , α 3 ), the Pauli matrices σ := (σ 1 , σ 2 , σ 3 ), and β = diag(1, 1, −1, −1). The operator on the right side of (4.23) is called the free Dirac operator. It acts on states ψ t (x) ∈ C 4 , called Dirac spinors. The underlying Hilbert space is L 2 (R 3 : C 4 ). The domain on which the free Dirac operator can be realized as a self-adjoint operator is H 1 (R 3 : C 4 ). The Foldy-Wouthuysen transform U FW allows to perform a block diagonalization, whereby the free Dirac operator takes the form This shows that the spectrum equals (−∞, −c 2 ] ∪ [c 2 , ∞). Physically, this means that states can possess "negative energy" and that there is an infinitely deep energy reservoir, the so-called Dirac sea. By adding electromagnetic fields and the charge conjugation operator, one can interpret states with negative energy as "antiparticles", i.e., particles with same mass but opposite charge. Such particles are called positrons.
Self-adjoint extensions of the Coulomb-Dirac operator. The one-particle Dirac operator describing the hydrogen atom can initially be defined on S(R 3 : C 4 ) and is formally given by the differential operator Scaling x → x/c and writing γ := Z/c shows that D H c,Z is unitarily equivalent to Weidmann [244] showed that D H γ is essentially self-adjoint on S(R 3 \ {0} : C 4 ) if and only if |γ| < √ 3/2, see also [231,Theorem 4.4]. For γ ∈ [ √ 3/2, 1] there is a "distinguished" (sometimes called "physically relevant") self-adjoint extension of D H γ . For γ ∈ ( √ 3/2, 1) it was established by Schmincke [197], Wüst [247] (see also Kalf, Schmincke, Walter, and Wüst [128] for a review of these results), Nenciu [181], and Klaus and Wüst [133]. According to Schmincke and Wüst, this realization stands out by the property that all states in the domain of the Coulomb-Dirac operator had finite potential energy. On the other hand, Nenciu's realization was distinguished by the fact that states had finite kinetic energy. Klaus and Wüst showed that both realizations coincide, and that the essential spectrum is (−∞, 1]∪ [1, ∞), see [134] (or [231, p. 117] for a textbook treatment). In summary, the domain dom(D H γ ) of the distinguished realization satisfies H 1 (R 3 : C 4 ), and the quadratic form domain is H 1/2 (R 3 : C 4 ). In particular, the expectation values of both kinetic and potential energy are finite in dom(D H γ ); this motivates the term "physically relevant extension". With the help of the sharp Hardy-Dirac inequality [42], Esteban and Loss [57] constructed a distinguished self-adjoint extension for γ = 1.
States in the domain of this operator need not have finite kinetic and potential energy separately. We remark that similar results in two dimensions (where the critical coupling is γ = 1/2) were proved by Warmt [242]. In this review we will only focus on the three-dimensional case and γ < 1.
Hydrogen eigenvectors and density. Similarly, the eigenvectors ψ D n,κ,m of the eigenvalue equation D H γ ψ D n,κ,m = λ n,κ ψ D n,κ,m are also well known, see, e.g., Pidduck [185] (in terms of Laguerre polynomials) or Gordon [96] and Darwin [32] (in terms of hypergeometric confluent functions). For nice pictorial representations, see, e.g., White [245] and, for textbook references, see Bethe [17,Formula (9.37)] and Thaller [232, p. 427]. The three-dimensional density for a Bohr atom for a given γ ∈ (0, 1) in spin-orbit channel κ ∈Ż is 36) and the total, three-dimensional hydrogenic density is These quantities are indeed well defined, as the following theorem demonstrates. To state it, let Moreover, for any ε > 0 there are constants A γ,ε , A γ > 0 such that for all x ∈ R 3 \ {0} one has The proof of this theorem is similar to that of Theorem 4.1.
Instability. In the Chandrasekhar model we call an atom "unstable", if the coupling constant is so large that the operator is unbounded from below. Recall that for γ < 1 the Coulomb-Dirac operator D H γ has a "distinguished" self-adjoint extension with the property that the expectation values of both kinetic and potential energy are finite in dom(D H γ ). Instability for the Coulomb-Dirac operator refers to the fact that all self-adjoint extensions of the Coulomb-Dirac operator have infinitely many eigenfunctions with infinite expectation value of the potential energy. This situation occurs when γ > 1, i.e., in the case when the lowest eigenvalue of D H γ has hit the lower threshold of the continuous spectrum. See, e.g., Hogreve [111, Theorem 2.1.(iii)] and the references therein for details and Thaller [231, p. 218] for an overview.
Brown-Ravenhall operator. According to Dirac, the "vacuum", i.e., the situation in which no (negatively charged) electrons with positive energies are present, is described by a completely filled negative energy continuum of the Dirac operator (the so-called "Dirac sea"), whereas only solutions to the free Dirac equation with positive kinetic energy should be regarded as "physical electrons". Dirac proposed the following interpretation [40, p. 362]: "The most stable states for an electron (the states of lowest energy) are those with negative energy and very high velocity. All the electrons in the world will tend to fall into these states with emission of radiation. The Pauli exclusion principle, however, will come into play and prevent more than one electron going into any one state. Let us assume there are so many electrons in the world that all the most stable states are occupied, or, more accurately, that all the states of negative energy are occupied except perhaps a few of small velocity. Any electrons with positive energy will now have very little chance of jumping into negative-energy states and will therefore behave like electrons are observed to behave in the laboratory. We shall have an infinite number of electrons in negative-energy states, and indeed an infinite number per unit volume all over the world, but if their distribution is exactly uniform we should expect them to be completely unobservable. Only the small departures from exact uniformity, brought about by some of the negative-energy states being unoccupied, can we hope to observe." Shortly after Dirac's equation was formulated, Breit [19,20,21] derived a relativistic wave equation for helium using the quantum electrodynamics of Heisenberg and Pauli, which reads As is explained in [22, p. 552] "the superscripts (1) and (2)  The system can make real transitions to states where one electron has a large negative energy and the other electron is in the positiveenergy continuum; thus equation (4.40) has no stationary solutions if interpreted in this way." In the language of spectral theory, the spectrum of a two-particle Coulomb-Dirac operator occupies the whole real axis (already without electron-electron repulsion), and all eigenvalues are embedded. This phenomenon is sometimes called Brown-Ravenhall disease [229,186,192]. As is well known from the non-relativistic theory, these are likely to turn into resonances when the electron-electron repulsion is turned on, which leads to unphysical consequences, such as the instability of the atom. We conclude this discussion by mentioning that despite these physical deficits, Oelker [182] recently showed that the single-particle Dirac operator may be extended to a self-adjoint multi-particle operator.
To remedy the above serious defects, Brown and Ravenhall proposed to only allow states with positive energy with respect to the free Dirac operator, i.e., they restricted the Hilbert space of admissible states to h c,0 := Λ c,0 (L 2 (R 3 : C 4 )) := 1 (0,∞) (−icα · ∇ + c 2 β)(L 2 (R 3 : C 4 )). (4.41) The energy of an electron in the Brown-Ravenhall picture is then given by The work [59] showed that this energy is bounded from below, if and only if γ = Z/c ≤ γ B with γ B := 2 π/2 + 2/π , (4.43) which implies Z < 125. For such γ the energy form can be extended to a closed quadratic form in h c,0 with form domain h c,0 ∩S(R 3 : C 4 ). The resulting self-adjoint operator constructed according to Friedrichs is called Brown-Ravenhall operator and is denoted by B c,Z . By scaling x → x/c the Brown-Ravenhall operator is seen to be unitarily equivalent to c 2 B 1,γ with γ = Z/c. We write B γ := B 1,γ . Although we will not use it in this review, we record the following convenient representation of the Brown-Ravenhall operator as a self-adjoint operator in L 2 (R 3 : we define the C 2×2 -valued functions Then the map maps L 2 (R 3 : C 2 ) unitarily onto h c,0 , cf. [59]. Therefore, B γ in h c,0 is unitarily equivalent to the operator with the "twisted potential" for any V : R 3 → C 4 , whenever meaningfully defined. From a technical point of view the representation (4.47) is important, e.g., in [59,23,91,168]. (Figuratively speaking, the transformation T c mollifies the Coulomb singularity and ensures the lower boundedness of B γ for coupling constants greater than 2/π. This transpires, e.g., in [91,Lemma 2.7] and [82,Lemma 5.2].) In [59] the authors showed that the Brown-Ravenhall operator B γ is bounded from below by −γ(π/4 − 1/π) − 1 when γ ≤ γ B . In fact, Tix [235,236] proved the stronger lower bound −γ B .
If γ < γ B , the essential spectrum of the Brown-Ravenhall operator B γ is [0, ∞) and the singular continuous spectrum is empty [59,Theorem 2]. Moreover, there are no embedded eigenvalues, and the spectrum in [0, ∞) is purely absolutely continuous [7,Theorem 3.4.1]. As in the Chandrasekhar case, the eigenvalues λ Z,n,j,ℓ of B c,Z in channel (j, ℓ) satisfy the bounds − Z 2 2(n + ℓ + 1) 2 ≥ λ Z,n,j,ℓ ≥ −const · where the constant in the second inequality can be chosen independently of γ. The lower bound is due to [91,Theorem 2.1], while the upper bound follows from the fact that the non-relativistic kinetic energy dominates the relativistic one. In particular, the Brown-Ravenhall eigenvalues are smaller than those of D H c,Z − c 2 due to the min-max-principle for operators with spectral gaps, cf. [99,98].
Furry operator. Naturally, the projection onto the positive spectral subspace of the free Dirac operator is not the only possibility to get rid of the positronic part of the wave functions. Furry and Oppenheimer [93] proposed rather to project onto the positive spectral subspace of the Coulomb-Dirac operator, i.e., the Hilbert space of admissible states is h c,Z := Λ c,Z (L 2 (R 3 : C 4 )) := 1 (0,∞) (D H c,Z )(L 2 (R 3 : C 4 )).
Since D H γ can be realized as a self-adjoint operator with form domain H 1/2 (R 3 : C 4 ) by Nenciu's method [181] when Z/c = γ < γ F with and dense in h c,Z . Therefore, the quadratic form is well-defined and bounded from below when γ ∈ (0, 1). By Friedrichs, there is a corresponding self-adjoint operator. The quadratic form domain of this operator is This operator is called Furry operator and denoted by F c,Z . Scaling x → x/c shows that the Furry operator is unitarily equivalent to c 2 F 1,γ with γ = Z/c and we write F γ := F 1,γ in the rescaled picture.
Mittleman operator. Although, at this point, the choice of the projections seems to be arbitrary and only justifiable by the comparison of the results with the measured quantities this is -according to Mittleman [171] not the case: the optimal projection and optimal ground state should be obtained by a minimax-principle, namely the infimum over the states in a class of fermionic Hilbert spaces defined by the positive spectral subspace of some Dirac operator −iα·∇+mc 2 β−ϕ followed by a supremum over a suitable class of potentials ϕ. However, this has not been implemented on a mathematical level. In fact there some elementary no-go results [9,11] that a potential mathematical implementation has to circumvent. For a more detailed review of Mittleman's principle and references, see, e.g., [54,Section 4.5].
For an associated ground state ψ of C N,V (or an approximate ground state on the Thomas-Fermi or Scott scale), we define the associated three-dimensional oneparticle ground state densities (4.53) Similarly as in the Schrödinger case, we define the one-dimensional angular momentum resolved version of ρ C for any ℓ ∈ N 0 by for r > 0. Note that ̺ ℓ,C (r), r > 0.  We record that Morozov and Vugalter [174] (see also Morozov [172], Jakubaßa-Amundsen [123] for HVZ theorems), and Matte and Stockmeyer [166] proved that for N = Z, the Brown-Ravenhall and Furry ground state energies E B/F (Z) := E B/F (Z, Z) are eigenvalues.
For an associated ground state ψ of E c,Z,N (or an approximate ground state on the Thomas-Fermi or Scott scale) in either the Brown-Ravenhall or Furry picture, we define the associated one-particle ground state densities (4.59) As in the Chandrasekhar case, we define a spin-orbit resolved version of ρ B/F . More precisely, for given spin-orbit coupling κ ∈Ż (recall (4.29)-(4.30)), we define  [228]. These operators are popular among quantum chemists, as they provide decent numerical results which are in good accordance with experimentally measured data. For instance, the Scott correction in the Furry picture (Formula (4.78) in Theorem 4.11) coincides astonishingly well with experimental data (see, e.g., [135]), see [100,Section 6] and [186,192] for textbook treatments.
(2) The Scott correction is also believed to be true when the mean field in the sense of Mittleman [171] is taken into account. However, this is so far only known in the Hartree-Fock approximation when the involved projection is given by the Dirac-Fock operator, see Fournais et al. [79].
In the following subsections we summarize results concerning the asymptotic expansion of the ground state energies and convergence of the one-particle ground state densities for the above-introduced models in the atomic and molecular cases. In particular we review the elements of the proof of the relativistic strong Scott conjecture for Chandrasekhar atoms in [87,85].

4.3.
Atoms without magnetic fields. As far as we know, the results below have only been proved in the neutral case N = Z. We believe that they also hold for ions.

Thomas-Fermi scale.
In all of the above relativistic models the leading order of the asymptotic expansion of the ground state energy is non-relativistic. These results indicate that electrons whose distances to the nucleus are of order Z −1/3 still behave non-relativistically, although they are being sucked into the nucleus. The following theorem was proved by Sørensen [184] (Chandrasekhar), by [23] and [91] Remarks 4.7 (Elements in the proof of Theorem 4.6). We make some remarks on Sørensen's proof [184] for the Chandrasekhar case. Since p 2 + 1 − 1 ≤ p 2 /2, it suffices to prove the lower bound, whose proof is similar to Lieb's simplified proof of Theorem 3.1, see [144,Theorem 5.1]. It can be split into the following steps.
(1) Reduce the linear many-particle problem to estimating the non-linear oneparticle quadratic form with the Thomas-Fermi density ρ TF Z and orthonormal orbitals {m ν } N ν=1 from below with the help of a correlation inequality (e.g., by Lieb and Oxford [183,142,149] or that of [160]).  [34], by [83] and the inequality p 2 + 1 − 1 ≥ |p| − 1. (Inequality (4.66) is often called Hardy-Lieb-Thirring inequality because of the homogeneity of the "unperturbed Hardy operator" |p| − 2/π |x| .) Here, V is a bounded function, supported on {|x| Z −o }. Using (4.66) one computes the contribution to the energy to be O(Z 2+3(1−o) ), which is more than Z 7/3 . For this reason another localization on the length scale Z −i is necessary. The Hardy-Lieb-Thirring inequality for V being supported on {|x| Z −i } then produces an o(Z 7/3 ) error when i ∈ (8/9, 1). The additionally introduced localization error can be controlled by an ε of kinetic energy, too.  Here semiclassical analysis is used. Roughly speaking, one compares the quantum energy to the classically expected energy The latter leads to the non-relativistic TF energy, since To that end a phase-space localization using coherent states [144, Theorem 5.1] is used. In fact, merely the localization errors coming from the phase-space localization force the position localization to the scale Z −1/3 . Theorem 4.6 is accompanied by the following convergence results for the ground state densities.
In all three formulae the convergence holds in Coulomb norm (see (2.5)-(2.6)) with convergence rate O(Z −3/16 ). In the Chandrasekhar case the convergence also holds when both sides are integrated against any U ∈ L 5/2 ∩ L 4 (R 3 ) and in the Brown-Ravenhall case when in addition U ∈ | · | −1 L ∞ is Lipschitz.

Scott scale.
As explained in the introduction, electrons in proximity of the nucleus are expected to generate relativistic effects that should be visible in the ground state energy and density on the spatial scale Z −1 . In fact, the Scott correction is relativistically lowered. The precise amount depends on the sum of the differences of the non-relativistic and the relativistic hydrogen eigenvalues. To that end, let λ S n , λ C n , λ B n , and λ F n denote the γ-dependent eigenvalues of the nonrelativistic operator S H γ in (1.9), of the Chandrasekhar operator C H γ in (4.2) (with q = 2), of the Brown-Ravenhall operator B γ . and of the Furry operator F γ . (The Furry eigenvalues coincide of course with those (4.35) of D H γ − 1.) We introduce the spectral shifts and record the following observation.  [90, p. 552], where it is shown that s C is monotone decreasing and finite. For s F the claim can be inferred from the explicitly known eigenvalues λ S n and λ F n , respectively. The continuity and monotonicity of s B follow from the explicit knowledge of the Schrödinger and Coulomb-Dirac eigenvalues, and the inequality λ F n ≥ λ B n . The following theorem concerning the energy asymptotics of Chandrasekhar atoms was proved independently by Solovej et al. [223], and the work [90] using different techniques. The result for Brown-Ravenhall atoms was proved by [91], and that for Furry atoms by [100].
In all of the above limits the error term can be quantified and is O(Z 47/24 ).
Remark 4.12. It is believed that these results also hold for ions (at least as long as the ionization degree is sufficiently small), since the electrons on length scales O(Z 0 ) should not disturb the energy generated by electrons on length scales O(Z −1 ). However, a rigorous proof is lacking.
The energetic results on the Scott scale are accompanied by recent results [87,85,170] for the density in the Chandrasekhar and Furry cases.
be arbitrary. Then the following statements hold.
(1) For the sake of clarity we restricted attention to the above class D of test functions, although the results actually hold for a substantially larger class.
(2) However, the exemplary test function class (4.79) is believed to be optimal. (a) Due to Kato's inequality, we cannot expect (at least not for ℓ = 0) (4.80)- (4.83) to hold for test functions, whose singularity is worse than |x| −1 .

4.3.3.
Elements in the proof of Theorem 4.13. The rest of this subsection is concerned with explaining the key elements of the proof of Theorem 4.13 in the Chandrasekhar case, i.e., the limits (4.80) and (4.82). For simplicity we set q = 1 here. We begin with the argument to prove (4.80) for a fixed angular momentum channel.
We follow the lines of Lieb and Simon [152], Baumgartner [12], and [117] by employing a linear-response argument. Let |ψ ψ| be a ground state density matrix of the atomic many-particle operator C Z (see (4.51)), and define, for U ∈ D and in slight abuse of notation, the perturbed operator Here Π ℓ is the orthogonal projection in L 2 (R 3 ) onto the ℓ-th angular momentum channel defined by and Π ℓ,ν acts as Π ℓ with respect to the ν-th particle. Since the singularity of U is Coulombic, C Z,λ is realized as a self-adjoint operator by Kato's inequality if γ < 2/π and |λ| is sufficiently small. Moreover, since (4.80) is linear in U we may assume U ≥ 0 and λ > 0 without loss of generality. By the linear response argument, the Scott correction (Theorem 4.11), and scaling with C H ℓ,γ as in (4.6). To compute the right side of (4.85) we have two options. (1) Find a majorant to apply the dominated convergence theorem to interchange lim inf λց0 and Tr. Then apply standard perturbation theory, i.e., the classical Hellmann-Feynman theorem for a single eigenvalue. This lead to the shorter proof in [85]. (2) Compute the derivative with respect to λ directly. This lead to the longer, original proof in [87]. Here we shall present the arguments of the longer proof, as we believe that it better unearths the involved mathematics of the relativistic strong Scott conjecture. Besides, it allows us to popularize a generalization of the classical Hellmann-Feynman theorem that may be of independent interest in the analysis of manyparticle problems.
First we state this generalized Hellmann-Feynman theorem with "natural" assumptions on the perturbation. However, this version is not applicable to our problem. Afterwards we state a generalization with weaker assumptions, which suffices for our purposes. Recall that an operator B is called relatively form trace class with respect to a self-adjoint operator A that is bounded from below, if (A + M ) −1/2 B(A + M ) −1/2 is trace class for some (and hence any) large enough M > 0. (4.86) In particular, S is differentiable at λ = 0, if and only if B| ker A = 0.
Remarks 4.16. (1) The relative form trace class assumption implies that the expression on the right of (4. 86), and consequently also that on the left, is finite.
(2) By the variational principle it follows that S is convex. Thus, S has left and right sided derivatives (cf. [217,Theorem 1.26]).
(3) If inf σ ess (A) > 0 or A − λB has only finitely many negative eigenvalues, then the derivative and the trace can be interchanged and the result follows from the classical Hellman-Feynman theorem. The point is that the formulae remain valid even when the bottom of the essential spectrum is zero, so that perturbation theory is not directly applicable.
In our application, A is the Chandrasekhar operator C H ℓ,γ (which has no zero eigenvalue) and B = U in L 2 (R + , dr). Thus, (4.80) would follow from Theorem 4.15, if one could verify its assumptions. By Kato's inequality (using γ < 2/π), one can replace A by the kinetic energy C ℓ . In this case, the relative trace class condition can be formulated explicitly using the Fourier-Bessel transform. But since (k + 1) −1 / ∈ L 1 (R + , dk), this shows that C H ℓ,γ cannot satisfy the relative trace class assumption, no matter how nice the perturbation B = U is.
For this reason, we proved and used a generalization of Theorem 4.15, where the relative trace class assumption is stated with respect to (A + M ) 2s for some s > 1/2. We will state this generalization in Theorem 4.18 in a moment. Using the following special case of [86, Theorem 1.1], the assumption of Theorem 4.18 can be recast as an assumption involving C ℓ instead of C H ℓ,γ . . Let 2/π ≥ γ ≥ 0, 0 < s < 3/2 − σ γ (with σ γ as in (4.16)), and V ∈ L 1 loc (R 3 ) satisfying −γ/|x| ≤ V ≤ 0. Then we have the quadratic form inequality |p| 2s (|p| + V ) 2s . (4.87) We are now ready to state the generalization of Theorem 4.15.  Thus, to conclude the proof of (4.80), we are left with showing the assumptions of Theorem 4.18. Since the test functions U in the above formulation of Theorem 4.13 are bounded by a multiple of the Coulomb potential, the condition (4.89) can be verified easily using Theorem 4.17. We now verify (4.88) with A = C H ℓ,γ and B = U .
Letting · 2 denote the Hilbert-Schmidt norm, we have, for any 1/2 < s ≤ 1, Clearly, functions U ∈ D, the test function space (4.79), satisfy U K (0) s < ∞. This concludes the sketch of the arguments in the proof of (4.80).
To prove the convergence of the total density c −3 ρ C (x/c), we also need to interchange the ℓ-summation with the limits Z → ∞ and λ → 0. This is done using the dominated convergence theorem. To apply it, it suffices to prove that there is an ε > 0 such that whenever 0 ≤ V ≤ γ/r. In our application, V is closely related to the Thomas-Fermi potential. For the proof of (4.90), see [87, Section 5].

4.4.
Molecules without magnetic fields. We immediately present the energetic result for Chandrasekhar molecules on the Scott scale, which was proved by Solovej et al. [223].
Theorem 4.21. Let the notations and assumptions be as in Theorem 4.19. In particular, fix K, z ∈ (0, 1) K , and r ∈ R 3K . Assume furthermore that there is γ 0 < 2/π such that max κ Z κ /c < γ 0 . Then  (1) Since c −1 |Z| is bounded, the prefactor c 2 /(8π 2 ) of the magnetic energy in the relativistic case is of order |Z| 2 (at least if we assume additionally that |Z|/c is kept constant when |Z|, c → ∞). This is much larger than in the non-relativistic case (Theorem 3.18), where c −2 |Z| was bounded. Thus, the selfgenerated magnetic field has to be much smaller in the relativistic case, which explains why it does not alter the Scott coefficient (contrary to the non-relativistic situation). In fact, the prefactor 8π in front of the field energy is irrelevant and can be replaced by any other fixed constant in the relativistic situation.
(2) The convergence of the density on the Scott scale has not been worked out so far.

Open questions
We collect some questions that -at least from our perspective -are interesting both from physical and mathematical points of view. Even more, one might ask whether the strong Scott conjecture would hold around each of the nuclei if (5.2) is true. A step in this direction was taken by Iantchenko et al. [117,Theorem 3] under the additional hypothesis that the minimal nuclear distance is bounded from below, e.g., by const |Z| − 1 4 . The analogue of E S can be defined for other many-particle models discussed in this review, e.g., those in Section 4. As far as we know, the Scott conjecture is also open for these problems. Again, we expect additivity of the energy up to order o(|Z| 5/3 ).
One could also consider a variant of this problem, where the kinetic energy of the nuclei is taken into account.
(2) It is folklore in quantum chemistry that chemical accuracy is achieved without taking a self-generated magnetic field into account. Therefore, we ask whether in physical models at least the Scott conjecture does not depend on the self-generated magnetic field. (3) The Scott conjectures for the energy and the density are open for no-pair operators where a self-consistent mean field is taken into account. In this case the Hilbert space of admissible one-particle states (in the x → x/c rescaled picture) is Λ χ (L 2 (R 3 : C 4 )) with Λ χ := 1 (0,∞) (D H γ + χ) with a (not necessary local) mean field χ. The resulting no-pair operator is sometimes called Fuzzy operator. Possible choices for χ are the Thomas-Fermi potential (the right sides of (2.11)-(2.12)), or the Hartree-Fock potential i |ϕ i | 2 * 1 | · | (x) − ϕ i (x)ϕ i (·) |x − ·| generated by a set of appropriately chosen orbitals {ϕ i } i . The latter choice is especially popular in quantum chemistry, see, e.g., [125,126,194,127]. It turns out, though, that numerically computed values of the ground state energy for Coulomb systems with Z 90 in the Fuzzy picture are quite close to those in the Furry picture; see also [171,192,195,100].
Nevertheless, from a mathematical point of view it is interesting to investigate the precise value of the Scott correction for the Fuzzy operator. The following tasks seem natural. (a) Show that for any (reasonable?) choice of the mean field χ, the leading order of the ground state energy is still the Thomas-Fermi energy. (b) Show the Scott conjecture for the Fuzzy model defined in the spirit of Mittleman (minimization of the electronic degrees followed by a maximization of the splitting). (c) Show that, to within the order of accuracy of Scott's correction, the maximization is attained in the Furry picture. (Recall that the Z 2correction in the Furry picture is exclusively generated by the effective one-particle problem involving the hydrogenic operator D H γ . Thus, our intuition is supported by the variational principle for operators with spectral gaps [99,173], which leads to the largest eigenvalues of Λ χ D H γ Λ χ , when one chooses Λ χ to be the Furry projection. This is also the underlying spirit of Schwinger's derivation [198] of the relativistic Scott correction.) (4) Can one show first and second order asymptotics for the Lieb-Loss model [140]? See Bach and Hach [6] for the first order in the non-relativistic setting. (5) Let us formulate some questions regarding the hydrogenic densities ρ H # with # ∈ {C, B, F }. For # = F one might be able to exploit the explicitly known eigenfunctions of the Coulomb-Dirac operator (cf. [96,32,185,17,231]) as in Theorem 3.11 by Heilmann and Lieb. as in Theorem 3.11? (Here ρ TF z=1 is the hydrogenic TF density with q = 2.) In case # = F one might be able to derive an asymptotic expansion as in (3.30) using the explicitly known eigenfunctions. (c) The non-relativistic density ρ H S (x) decreases monotonically in |x| > 0. Is this also true for the relativistic hydrogenic densities? (6) As we have seen in Theorem 2.10 and the ensuing discussion, Weizsäcker's parameter A in the TFW functional (2.23) can be tuned to achieve either energy agreement, or agreement of the TFW density near the origin. This leads to the following (vague) question of whether one can construct a modification of Weizsäcker's gradient term that gives simultaneously the energy and the density at the origin correctly.