How to place optimally four electrons on a sphere? $\newcommand{\S}{\mathbb{S}^2}$
Let $x_1,x_2,x_3,x_4 \in \mathbb{S}^2$ be points on the unit sphere, that minimizes the quantity
$$
E(x_1,x_2,x_3,x_4)=\sum_{i < j}\frac{1}{\| x_i - x_j \|},
$$
where $\| x_i - x_j \|$ denotes the Euclidean distance in $\mathbb{R}^3$.
$E$ is the electrostatic potential energy of $4$ electrons constrained to lie on the unit sphere.


It is claimed in various sources that the $x_i$ are the vertices of a regular tetrahedron. (see e.g. Wikipedia).

Question: I would like to find a reference for a proof of that fact. (Or a self-contained proof produced here, if it is not too large). I found various sources on this problem, known as the Thomson problem, but could not actually find a paper containing a proof of the four-particle case.

Edit:
Here is a proof that the tetrahedron $(x_1,x_2,x_3,x_4)$ is at a local minimum:
Let $\beta(t)$ be a path in $\mathbb{S}^2$, $\beta(0)=x_4, \dot \beta(0)=w \in T_{x_4}\S$.
Consider the path $\alpha(t)=(x_1,x_2,x_3,\beta(t))$. Since $\alpha$ changes only the fourth component, we get
$$
\frac{d}{dt}E(\alpha(t))=\sum_{i=1}^3\frac{d}{dt}| x_i - \beta(t) |^{-1}.
$$
$$
\frac{d}{dt}| x_i - \beta(t) |^{-1}=\frac{d}{dt}(| x_i - \beta(t) |^2)^{-\frac{1}{2}}=| x_i - \beta(t) |^{-3}\langle x_i,\dot \beta(t)\rangle. \tag{1}
$$
Differentiating again, we get
$$
\frac{d^2}{dt^2}| x_i - \beta(t) |^{-1}=| x_i - \beta(t) |^{-3}\big(3| x_i - \beta(t) |^{-2}\langle x_i,\dot \beta(t)\rangle^2+\langle x_i,\ddot \beta(t)\rangle \big)\tag{2}.
$$
In particular, for a tetrahedron, we get
$$
dE_p(0,0,0,w)=\sum_{i=1}^3 | x_i - x_4 |^{-3}\langle x_i,w\rangle=a^{-3}  \langle \sum_{i=1}^3 x_i,w\rangle= 
\langle -x_4,w \rangle=0,
$$
where we used the facts that $\sum_{i=1}^4 x_i=0$, and $\langle x_4,w \rangle$, since $x \in T_{x_4}\S$.
Here $a=\sqrt\frac{8}{3}$ is the edge length of the tetrahedron. Moreover,
$$
a^{3}\frac{d^2}{dt^2}|_{t=0}E(\alpha(t))=3a^{-2}\sum_{i=1}^3\langle x_i,w\rangle^2+\langle \sum_{i=1}^3x_i,\ddot \beta(0)\rangle.
$$
Differentiating twice $\langle \beta(t),\beta(t) \rangle=1$, we get $ \langle x_4,\ddot \beta(0) \rangle=-|w|^2$, so
$$
\frac{d^2}{dt^2}|_{t=0}E(\alpha(t))=3a^{-2}\sum_{i=1}^3\langle x_i,w\rangle^2+|w|^2 \ge 0,
$$
and is strictly positive whenever $w \neq 0$. This proves that the tetrahedron is indeed a point of local minimum.
 A: Let

*

*$f : (0,4] \to (0,\infty)$ be any monotonic decreasing convex function.

*$p_1,\ldots, p_4$ be any $4$ distinct points on $S^2$.

*$q_1,\ldots, q_4$ be $4$ points on $S^2$ forming a regular tetrahedron.

Notice
$$\begin{align}\sum_{i<j} |p_i - p_j|^2 
&= \frac12\sum_{i,j}|p_i - p_j|^2
= \frac12\left(4\sum_j |p_j|^2 + 4\sum_i|p_i|^2 - 2\sum_{ij} p_i\cdot p_j\right)\\
&= 16 - \left|\sum_i p_i\right|^2
\end{align}$$
Since $f$ is convex.
$$\sum_{i<j}f\left(|p_i - p_j|^2\right) \ge 6f\left(\frac16\sum_{i<j}|p_i-p_j|^2\right) = 6f\left[\frac83-\frac16\left|\sum_ip_i\right|^2\right]$$
Since $f$ is monotonic decreasing,
$$\sum_{i<j}f\left(|p_i - p_j|^2\right) \ge 6f\left(\frac83\right)$$
Notice $\displaystyle\;|q_i - q_j|^2 = \frac83\;$ for all $i \ne j$, we have
$$\sum_{i<j}f\left(|p_i - p_j|^2\right) \ge \sum_{i<j}f\left(|q_i - q_j|^2\right)$$
Now the map $\displaystyle\;x \mapsto \frac{1}{\sqrt{x}}\;$ is monotonic decreasing and convex on $(0,4]$. This leads to
$$\sum_{i<j} \frac{1}{|p_i - p_j|} \ge \sum_{i<j} \frac{1}{|q_i - q_j|}$$
As a result, the electrostatic potential energy is minimized when $p_1,\ldots,p_4$ are vertices of a regular tetrahedron.
A: Nice question! I'm sure there are better ways but I just came up with this approach.
Define the origin as the centre of the sphere and construct a regular tetrahedron on the sphere with its vertices at  $\mathbf{a} , \mathbf{b} , \mathbf{c} , \mathbf{d}$. We will show that $E$ is a local minimum with the charges at these vertices.
Note that, by symmetry, $\mathbf{a} + \mathbf{b} + \mathbf{c} + \mathbf{d}= \mathbf{0} $.
Suppose we move the charge at $\mathbf{d}$ by a small vector $\mathbf{\delta}$. We will calculate $\delta E$, the change in $E$ when we make this change.
$$E=\|\mathbf{a} - \mathbf{b} \|^{-1} + \|\mathbf{a} - \mathbf{c} \|^{-1} + \|\mathbf{a} - \mathbf{d} \|^{-1} + \|\mathbf{b} - \mathbf{c} \|^{-1} + \|\mathbf{b} - \mathbf{d} \|^{-1} + \|\mathbf{c} - \mathbf{d} \|^{-1} \tag1$$
$$\delta E=\sum_{ \mathbf{r}=\mathbf{a},\mathbf{b},\mathbf{c} }  \|\mathbf{r} - \mathbf{d}  - \mathbf{\delta}  \|^{-1} - \|\mathbf{r} - \mathbf{d} \|^{-1}  \tag2$$
But we can see that
$$\|\mathbf{r} - \mathbf{d} \|^{-1} =(\mathbf{r}^2 + \mathbf{d}^2  -2 \mathbf{r} \cdot \mathbf{d} )^{-\frac12} \tag3 $$
and
$$\|\mathbf{r} - \mathbf{d} - \mathbf{\delta}   \|^{-1}  =(\mathbf{r}^2 + \mathbf{d}^2 +\mathbf{\delta}^2  -2 \mathbf{r} \cdot \mathbf{d}  -2 \mathbf{r} \cdot \mathbf{\delta}  +2 \mathbf{d} \cdot \mathbf{\delta}  )^{-\frac12} \tag 4 $$
But $\mathbf{d} \cdot \mathbf{\delta}$ must be zero (because we are restricted to moving the charge on the surface of the sphere and the surface is perpendicular to the radius from the origin to point $\mathbf{d}$) and $\mathbf{\delta}^2$ can be ignored, so
$$\|\mathbf{r} - \mathbf{d} - \mathbf{\delta}   \|^{-1} =(\mathbf{r}^2 + \mathbf{d}^2   -2 \mathbf{r} \cdot \mathbf{d}  -2 \mathbf{r} \cdot \mathbf{\delta}   )^{-\frac12} $$
$$=(\mathbf{r}^2 + \mathbf{d}^2   -2 \mathbf{r} \cdot \mathbf{d})^{-\frac12}\left(1  -2 \frac{\mathbf{r} \cdot \mathbf{\delta} }{\mathbf{r}^2 + \mathbf{d}^2   -2 \mathbf{r} \cdot \mathbf{d}}  \right)^{-\frac12} $$
(ignoring terms in $\delta^2$):
$$=(\mathbf{r}^2 + \mathbf{d}^2   -2 \mathbf{r} \cdot \mathbf{d})^{-\frac12}\left(1  + \frac{\mathbf{r} \cdot \mathbf{\delta} }{\mathbf{r}^2 + \mathbf{d}^2   -2 \mathbf{r} \cdot \mathbf{d}}  \right) \tag 5 $$
and so, from (2), (3) and (5):
$$\delta E=\sum_{ \mathbf{r}=\mathbf{a},\mathbf{b},\mathbf{c} }  
(\mathbf{r}^2 + \mathbf{d}^2   -2 \mathbf{r} \cdot \mathbf{d})^{-\frac32}(\mathbf{r} \cdot \mathbf{\delta}) $$
$$\delta E=\sum_{ \mathbf{r}=\mathbf{a},\mathbf{b},\mathbf{c} }  
\|\mathbf{r}  + \mathbf{d}\|^{-3}(\mathbf{r} \cdot \mathbf{\delta})  \tag6$$
But it must be that $\|\mathbf{a}  + \mathbf{d}\|=\|\mathbf{b}  + \mathbf{d}\|=\|\mathbf{c}  + \mathbf{d}\|=k$ by the symmetry of the regular tetrahedron. So
$$\delta E=k^{-3}\sum_{ \mathbf{r}=\mathbf{a},\mathbf{b},\mathbf{c} }  
(\mathbf{r} \cdot \mathbf{\delta})= k^{-3}(\mathbf{a}+\mathbf{b}+\mathbf{c}) \cdot \mathbf{\delta}= k^{-3}(-\mathbf{d}) \cdot \mathbf{\delta} =\mathbf{0}\tag 7$$
And so $\delta E=0 $ and so the regular tetrahedral layout is a local minimum. I suppose it could actually be a maximum but I'm going to leave it there!
