---
abstract: 'We determine the ground state energy of atoms and quantum dots whose number $N$ of electrons is large. We show that the dominant terms of the energy are those given by a semiclassical Hartree-Fock theory. Correlation effects appear at the order $N\ln N$ for atoms and the order $N$ for quantum dots. We compute them. The semiclassical Hartree-Fock theory creates oscillations in the ground state energy as a function of $N$. These oscillations reflect the dynamics of a classical particle moving in the presence of the Thomas-Fermi potential. The dynamics is regular for atoms and some dots, but we present the case of a dot where this dynamics is fully chaotic and we compute the oscillating part of the ground state energy in this case.'
author:
- Hervé Kunz
- Rico Rueedi
title:  Ground state energy of large atoms and quantum dots 
---

Since the introduction by Thomas [@Thomas-PCPS23] and Fermi [@Fermi-RANL6] of the statistical theory of atoms, many efforts have been made to systematically improve on it. For neutral atoms with $N$ electrons, an expansion of the ground state energy, whose first term is Thomas-Fermi (TF) theory, was determined up to the order $N^{5/3}$ [@Morgan-1996]. This expansion was based on Hartree-Fock (HF) theory. But after Lieb and Simon [@LiebSimon-AdvMath23] proved that TF theory becomes exact when $N\to\infty$, part of the expansion was put on a rigorous basis [@FeffermanSeco-HPA71]. A crucial step further, inside HF theory, was made by Schwinger and Englert [@EnglertSchwinger-PRA32b] who showed that there exist oscillating corrections to the ground state energy, of order $N^{4/3}$, and period $N^{1/3}$. They interpreted such terms as indicating shell effects.

We consider quantum dots as artificial atoms made of $N$ electrons, confined by a potential and constrained to move in two dimensions. The determination of the dependence of the ground state energy of quantum dots as a function of $N$ is of great interest. It was proven by Lieb, Solovej, and Yngvason [@LiebSolovejYngvason-PRB51] that in this case also, TF theory becomes exact when $N\to\infty$. This suggests that like in the case of atoms, corrections to it should be determined. An important difference with atoms is the fact that the confining potential isn’t necessarily rotationally symmetric, and can be changed.

We have found a general formula for the ground state energy of both atoms and dots, appropriate to the case when $N{^{-1}}$ is small. It is based on the fact that in this case, we are confronted to a semiclassical mean field theory corrected by fluctuation effects. In this way, we have shown that HF theory is valid, up to a certain order in $N{^{-1}}$ where correlation effects become important. We have determined the dominant correlation correction, which in the case of atoms, is $N\ln N$.

We show that the oscillating contributions to the energy appearing in HF theory can be interpreted as coming from periodic orbits of a classical particle moving in the TF potential. In atoms, this motion is integrable. But in dots, the confining potential can be chosen so that the motion is fully chaotic, and we compute the ground state energy for a specific example of such a chaotic quantum dot.

We can put non-relativistic atoms and quantum dots on the same footing, if we describe them, in appropriate units, by the $d$-dimensional Hamiltonian of $N$ electrons $${\bar{H}}=
-\sum_{j=1}^N{\bar{{\Delta}}}_j
+Z\sum_{j=1}^N{V}({\bar{{{\bm{x}}}}}_j)
+{\frac{1}{2}}\sum_{j\neq k}^N
\frac{1}{{|{\bar{{{\bm{x}}}}}_j-{\bar{{{\bm{x}}}}}_k|}}.
\label{def:H:srv}$$ For atoms ($d=3$), $Z$ is the number of protons and ${V}({\bar{{{\bm{x}}}}})=-\smash{{|{\bar{{{\bm{x}}}}}|}{^{-1}}}$. For dots ($d=2$), $Z$ will be identified with $N$, a possible dependence of ${V}({\bar{{{\bm{x}}}}})$ on $N$ being neglected at this stage.

For the purpose of considering systems where $Z$ (and $N$) is large, we rescale the coordinates so that ${\bar{H}}=Z^{2-2/d}H$, where $$H =
-{\epsilon}^2\sum_{j=1}^N{\Delta}_j
+\sum_{j=1}^N{V}({{\bm{x}}}_j)
+\frac{1}{2Z}\sum_{j\neq k}^N
\frac{1}{{|{{\bm{x}}}_j-{{\bm{x}}}_k|}},
\label{def:H}$$ the parameter ${\epsilon}$, playing the role of $\hbar$, being ${\epsilon}\doteq Z^{-1/d}$. In this form, the Coulomb interaction is of mean field type, and, when $Z$ is large, the Hamiltonian has to be considered in the semiclassical limit. Our strategy to determine the ground state energy $E$ of $H$ is then the following: considering first ${\epsilon}$ and $Z$ as independent, we compute an asymptotic expansion in $Z{^{-1}}$ of the pressure ${P}$ of the Hamiltonian $H$, then take the zero temperature limit $({\beta}\to\infty$), and finally look at the semiclassical limit ${\epsilon}\to0$.

With this strategy in mind we first make a Hubbard-Stratonovich transformation on the Coulomb interaction term in the coherent state path integral representation of the partition function ${Q}$ [@NegeleOrland-1988]. After an appropriate shift, the partition function can be written as the average value over an auxiliary Gaussian field $\theta$ $${Q}={\Bigl<e^{ZS(\theta)}\Bigr>_{\theta}},$$ $S(\theta)$ being an action such that its saddle point, which dominates ${Q}$ when $Z$ is large, corresponds to the Hartree approximation.

The final exact result for ${E}$ can be put in the form $${E}={E_\textsc{h}}+{E_\textsc{x}}+{E_{\textsc{c};1}}+{E_{\textsc{c};2}}.$$ The first term is the Hartree energy, given by $${E_\textsc{h}}\doteq
{\mu}N-\frac{Z}{2}\int
\frac{{d {{\bm{x}}}}{d {{\bm{y}}}}}{{|{{\bm{x}}}-{{\bm{y}}}|}}
{\rho_0}{({{\bm{x}}})}{\rho_0}{({{\bm{y}}})}-2\int^{\mu}{d e{\,}}{D}(e).$$ ${D}(e)=\operatorname{tr}{{\vartheta}(e-h)}$, with ${\vartheta}$ the Heaviside step function, is the integrated density of states of the one-body Hamiltonian $h\doteq-{\epsilon}^2{\Delta}+{W}{({{\bm{x}}})}$, with ${W}{({{\bm{x}}})}$ the self-consistent potential given by $${W}{({{\bm{x}}})}\doteq
{V}{({{\bm{x}}})}+
\int\frac{{d {{\bm{y}}}}}{{|{{\bm{x}}}-{{\bm{y}}}|}}{\rho_0}{({{\bm{y}}})}.$$ ${\rho_0}{({{\bm{x}}})}$ is the density ${\rho_0}{({{\bm{x}}})}\doteq 2{n}{({{\bm{x}}}|{{\bm{x}}})}/Z$ , where ${n}{{({{\bm{x}}}|{{\bm{y}}})}}\doteq {\langle{{\bm{x}}}|{{\vartheta}({\mu}-h)}|{{\bm{y}}}\rangle}$ is the density matrix associated to $h$. Finally, the chemical potential ${\mu}$ is fixed by the constraint $\partial_{\mu}{E}=0$. The second term is the exchange energy $${E_\textsc{x}}\doteq
-\frac{1}{Z}\int\frac{{d {{\bm{x}}}}{d {{\bm{y}}}}}{{|{{\bm{x}}}-{{\bm{y}}}|}}
{n}^2{{({{\bm{x}}}|{{\bm{y}}})}}.$$ This term also appears in HF theory, where it is the dominant exchange contribution when $Z$ is large [@Englert-1988]. The remaining terms describe correlation effects, with $${E_{\textsc{c};1}}\doteq
{\lim_{{\beta}\to\infty}}-\frac{1}{2{\beta}}\sum_{n=2}^{\infty}
\frac{1}{nZ^n}
\operatorname{Tr}(C\Gamma)^n,
\label{def:ECR1}$$ and $${E_{\textsc{c};2}}\doteq
{\lim_{{\beta}\to\infty}}-\frac{1}{{\beta}}
\ln
{\Bigl<\exp\Big[
-2\sum_{n=3}^{\infty}\frac{(-i)^n}{nZ^{n/2}}
\operatorname{Tr}(K\theta)^n\Big]\Bigr>_{\theta}}.$$ ${C}$ is the operator of kernel $${C}{({{\bm{x}}},{\tau_x}|{{{\bm{y}}}},{\tau_y})}=
\frac{1}{{|{{{\bm{x}}}}-{{{\bm{y}}}}|}}{{\delta}\big({\tau_x}-{\tau_y}\big)},$$ $0\leq{\tau_x},{\tau_y}\leq{\beta}$ being imaginary times. $K$ is the operator $$K=
[\partial_\tau-{\mu}+h]{^{-1}},$$ $\Gamma$ is the operator of kernel $$\Gamma{({{{\bm{x}}}},{\tau_x}|{{{\bm{y}}}},{\tau_y})}=
K{({{{\bm{x}}}},{\tau_x}|{{{\bm{y}}}},{\tau_y})}
K{({{{\bm{y}}}},{\tau_y}|{{{\bm{x}}}},{\tau_x})},$$ and finally, $\theta({{\bm{x}}},\tau)$ is a Gaussian field of covariance $[{C}{^{-1}}-\Gamma/Z\big]{({{{\bm{x}}}},{\tau_x}|{{{\bm{y}}}},{\tau_y})}$.

We recover TF theory by taking the semiclassical limit in the Hartree energy. Indeed, in this limit, the density ${\rho_0}{({{\bm{x}}})}$ and the integrated density of states $D(e)$ become $${\rho_0}{({{\bm{x}}})}=
\gamma_d{\mu}_+^{d/2}{({{\bm{x}}})},$$ where ${\mu}_+{({{\bm{x}}})}\doteq({\mu}-{W}{({{\bm{x}}})}){{\vartheta}({\mu}-{W}{({{\bm{x}}})})}$, and $${D}(e)=
\gamma_d{\int{d {{\bm{x}}}{\,}}}
\big(e-{W}{({{\bm{x}}})}\big)_+^{d/2},$$ where $\gamma_2=(2\pi){^{-1}}$ and $\gamma_3=(3\pi^2){^{-1}}$, and ${\mu}$ is fixed by the constraint $${\int{d {{\bm{x}}}{\,}}}{\rho_0}{({{\bm{x}}})}=1.$$

At the next stage, we compute the correlation energies. For ${E_{\textsc{c};1}}$, we simply need to take the limit ${\epsilon}\to0$. For ${E_{\textsc{c};2}}$, we first compute the terms up to the order $Z^{-3}$, then take the semiclassical limit. We recall that the true ground state energy is given by ${{\bar{E}}}=Z^{2-2/d}E$.

In the case of neutral atoms, and moving to standard atomic units, we find $${{\bar{E}}_{\textsc{c};1}}=
-0.062\ N\ln N-
[0.034+2x_2]N+
{\mathcal{O}(N^{2/3})},$$ where $$x_2\doteq
\frac{1}{2}
\sum_{n=3}^\infty
\frac{1}{n}
\int_{-\infty}^\infty\frac{{d \Omega}}{2\pi}\operatorname{tr}T_\Omega^n,$$ $T_\Omega$ being certain operator on $L^2({\mathbb{R}}^3)$, and $${{\bar{E}}_{\textsc{c};2}}=0.98\ N+{\mathcal{O}(N^{2/3})}.$$ It is remarkable that the dominant term for the correlation energy $-0.062\ N\ln N$ is the same as the logarithmically divergent second order correlation energy of the jellium [@Gell-MannBrueckner-PR106], an effective cut-off ${\epsilon}$ on some momentum making the integrals finite in the case of atoms. Other parts of the correlation energy are the same as in the jellium. For its correlation energy, the atoms looks like a jellium with a *local* Fermi surface defined by ${{\bm{p}}}^2\leq{\mu}_+{({{\bm{x}}})}$.

-\#1[\#1]{}(\#1,\#2)\#3[\[\]\[\]]{}\[cc\]\[cc\][$-\bar{E}_\textsc{cr}/N$ (Ha)]{}\[cl\]\[cl\]\[cl\]\[cl\]\[cc\]\[cc\][$N$]{}\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]![Per electron correlation energy of neutral atoms. Data points correspond to experimental values [@ChakravortyDavidson-JPC100] and numerical values obtained with an extended HF approach [@ClementiCorongiu-IJQC62], available for atoms containing up to $18$ and $55$ electrons, respectively. For sufficiently large atoms, the contribution ${{\bar{E}}_\textsc{c}}/N=-0.062\ln N-0.18$ deviates from experimental and extended HF values essentially by an oscillating term.[]{data-label="A:CR"}](f1.eps "fig:"){width="68mm"}

To summarize, for neutral atoms we have found the correlation contribution ${{\bar{E}}_\textsc{c}}\doteq{{\bar{E}}_{\textsc{c};1}}+{{\bar{E}}_{\textsc{c};2}}$ to the ground state energy $${{\bar{E}}_\textsc{c}}=-0.062\ N\ln N+cN+{\mathcal{O}(N^{2/3})},
\label{hECR}$$ and not, as has been suggested, a $N^{4/3}$ contribution [@Englert-1988; @ClementiCorongiu-IJQC62]. In figure \[A:CR\], a comparison is made between experimental [@ChakravortyDavidson-JPC100] and extended Hartree-Fock [@ClementiCorongiu-IJQC62] values for the correlation energy, and formula taking $c=-0.018$ for the best fit. We have not computed numerically $x_2$. In the case of quantum dots, and in the units of , we find $${{\bar{E}}_{\textsc{c};1}}=-0.15\ N+{\mathcal{O}(N^{1/2}\ln N)},$$ and $${{\bar{E}}_{\textsc{c};2}}=
-0.15\ N
+ 2c_1N
+{\mathcal{O}(N^{1/2})},$$ where $$c_1\doteq
\frac{1}{2\pi^4}
\int\frac{{d {{\bm{x}}}}{d {{\bm{y}}}}}{{|{{\bm{x}}}-{{\bm{y}}}|}}
{\mu}^{1/2}_+{({{\bm{x}}})}{\mu}^{1/2}_+{({{\bm{y}}})}.$$ It can be noticed that the dominant contribution to ${E_{\textsc{c};1}}$ is universal, *i.e.* independent of the confining potential.

Let us now look at the Hartree-Fock type contribution ${{\bar{E}}_\textsc{hf}}={{\bar{E}}_\textsc{h}}+{{\bar{E}}_\textsc{x}}$, which needs to be computed semiclassically. ${{\bar{E}}_\textsc{h}}$ depends essentially on the density of states ${d}(e)=\partial_e{D}(e)$. It is standard now to decompose ${d}(e)$ into two parts, semiclassically, as ${d}(e)={{d}_\text{s}}(e)+{{d}_\text{osc}}(e)$ [@BrackBhaduri-1997]; and, accordingly, to decompose ${{\bar{E}}_\textsc{hf}}$ as ${{\bar{E}}_\textsc{hf}}={{\bar{E}}_{\textsc{hf};\text{s}}}+{{\bar{E}}_{\textsc{hf};\text{osc}}}$. The smooth part ${{d}_\text{s}}(e)$ is given by an asymptotic expansion in ${\epsilon}$, whose coefficients are some integrals depending on ${W}{({{\bm{x}}})}$. But the oscillating part $${{d}_\text{osc}}(e)=
\sum_{\gamma}A_{\epsilon}(e,\gamma)\cos\Big(\frac{1}{{\epsilon}}S(e,\gamma)+\sigma_\gamma\frac{\pi}{2}\Big)$$ is given by a sum over the periodic orbits $\gamma$ of a classical particle moving in the potential ${W}{({{\bm{x}}})}$, where $S(e,\gamma)$ is the classical action along the orbit, $\sigma_\gamma$ is the orbit’s Maslov index, and $A_{\epsilon}(e,\gamma)$ depends on the orbit’s period and stability. ${{d}_\text{osc}}(e)$ therefore depends crucially on the nature of the dynamics associated to ${W}{({{\bm{x}}})}$.

In the case of atoms, the smooth part is known up to the order ${\epsilon}{^{-1}}$, and further corrections both to the Hartree and exchange energies would be required to reach the order ${\epsilon}$ of the correlation energy. The semiclassical determination of these corrections is particularly delicate due to the Coulomb singularity of ${W}{({{\bm{x}}})}$ at the origin, and the slow decay of the Coulomb potential appearing in the exchange term. Therefore, we have not undertaken this task. The oscillating part was computed by Schwinger and Englert, who did not give a dynamical interpretation of it. The spherically symmetric potential ${W}{({{\bm{x}}})}$ makes the dynamics integrable in this case.

In the case of dots, we have completely determined the smooth part of ${{\bar{E}}_\textsc{hf}}$ at the same order as the correlation energy. It is given by $$\begin{gathered}
{{\bar{E}}_{\textsc{hf};\text{s}}}=
{\mu}N^2
-\frac{N^2}{4\pi}
{\int{d {{\bm{x}}}{\,}}}{\mu}_+{({{\bm{x}}})}-\frac{N^2}{8\pi^2}
\int\frac{{d {{\bm{x}}}}{d {{\bm{y}}}}}{{|{{\bm{x}}}-{{\bm{y}}}|}}
{\mu}_+{({{\bm{x}}})}{\mu}_+{({{\bm{y}}})}\\
-\frac{2N^{3/2}}{3\pi^3}
{\int{d {{\bm{x}}}{\,}}}
{\mu}_+^{3/2}{({{\bm{x}}})}+\frac{N}{24\pi}
{\int{d {{\bm{x}}}{\,}}}{\Delta}{W}{({{\bm{x}}})}{{\vartheta}\big({\mu}{({{\bm{x}}})}\big)}+c_2N,\end{gathered}$$ where $$c_2\doteq
\frac{1}{\pi^3}
\frac{\displaystyle
\Big[
{\int{d {{\bm{x}}}{\,}}}{\mu}_+^{1/2}{({{\bm{x}}})}a{({{\bm{x}}})}\Big]^2
}{\displaystyle
{\int{d {{\bm{x}}}{\,}}}{{\vartheta}\big({\mu}{({{\bm{x}}})}\big)} a{({{\bm{x}}})}},$$ $a{({{\bm{x}}})}$ being the solution of the integral equation $$a{({{\bm{x}}})}=
1-\frac{1}{2\pi}
\int\frac{{d {{\bm{y}}}}}{{|{{\bm{x}}}-{{\bm{y}}}|}}
{{\vartheta}\big({\mu}{({{\bm{y}}})}\big)}a{({{\bm{y}}})};$$ and where $${W}{({{\bm{x}}})}=
{V}{({{\bm{x}}})}+
\frac{1}{2\pi}
\int\frac{{d {{\bm{y}}}}}{{|{{\bm{x}}}-{{\bm{y}}}|}}{\mu}_+{({{\bm{y}}})}.$$

-\#1[\#1]{}(\#1,\#2)\#3[\[\]\[\]]{}\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[tc\]\[cc\]\[cc\][$\bar{E}_{\textsc{osc}}/N^{1/2}$ ($10^{-2}$)]{}\[Bc\]\[Bc\][$N^{1/2}$]{}\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]\[cr\]![Oscillating contribution to the ground state energy for the chaotic quantum dot of quartic confinement ($\lambda=1.2$), computed with the 30 shortest periodic orbits of the classical system.[]{data-label="E:OS"}](f2.eps "fig:"){width="68mm"}

The oscillating corrections to the ground state energy are of order $N^\alpha$, and of period $N^{1/2}$, where $\alpha=3/4$ in the integrable case and $\alpha=1/2$ in the fully chaotic case. Note that in both cases and contrary to the case of atoms, the oscillating corrections are less important than correlation effects. If ${V}{({{\bm{x}}})}$ is a harmonic well, the dynamics induced by ${W}{({{\bm{x}}})}$ is integrable and the oscillating part has been computed [@Plocek-2006]. It is particularly interesting to consider a case where the dynamics induced by ${W}{({{\bm{x}}})}$ is fully chaotic. For this purpose we have chosen the quartic potential ${V}{({{\bm{x}}})}={\alpha}(x_1^4/{\beta}+x_2^4{\beta}-\lambda x_1^2x_2^2)$, with ${\alpha}=10^{-3}$ and ${\beta}=\pi/4$ [@BohigasTomsovicUllmo-PhysRep223; @BrackEA-JPhysA36]. The nature of the dynamics of ${W}{({{\bm{x}}})}$ depends on the coupling constant ${\lambda}$, and when $\lambda=1.2$, the dynamics is fully chaotic (whereas that of ${V}{({{\bm{x}}})}$ is not) [@Rueedi-2009]. In this case $${{\bar{E}}_{\textsc{hf};\text{s}}}=
 0.14\ N^2
-0.067\ N^{3/2}
+0.053\ N$$ and ${{\bar{E}}_{\textsc{hf};\text{osc}}}$, computed by summing over the 30 shortest periodic orbits, is shown in figure \[E:OS\].

To conclude, we have found a systematic way to compute the ground state energy of atoms and quantum dots when the number of electrons is large. A semiclassical HF theory has been shown to give the dominant terms. But correlation effects have been computed which determine the validity of HF theory. We have also shown that shell effects result from the dynamics of a classical particle moving in the presence of the TF potential. This allows, in the case of dots, to see the impact on the oscillations of the energy of a regular or chaotic motion of this classical particle.

We thank C. Plocek for useful and stimulating discussions at the early stage of this work. This work was supported by the Fonds National Suisse de la Recherche Scientifique.

[18]{} natexlab\#1[\#1]{}bibnamefont \#1[\#1]{}bibfnamefont \#1[\#1]{}citenamefont \#1[\#1]{}url \#1[`#1`]{}urlprefix\[2\][\#2]{} \[2\]\[\][[\#2](#2)]{}

, ****, ().

, ****, ().

, in **, edited by (, ), chap. , pp. .

, ****, (). , ****, ().

, ****, (). , , , ****, (). , **, Frontiers in Physics (, ).

, **, Lecture Notes in Physics (, ).

, ****, (). , ****, (). , ****, (). , **, Frontiers in Physics (, ).

, Ph.D. thesis, (), .

, , , ****, (). , , , , ****, (). , Ph.D. thesis, (), .
