Solution 1:

On the one hand,

$$\frac{1}{r^{n-1}}\frac{d}{dr}\left(r^{n-1}f'\right) = \frac{1}{r^{n-1}} \left( (n-1)r^{n-2} f' + r^{n-1} f'' \right) = \frac{n-1}{r}f' + f''$$

On the other, the gradient of $f$ is $\displaystyle \nabla f = \frac{\vec{r}}{r}f'$ and thus the Laplacian of $f$ is

$$\nabla^2 f = \left( \frac{n}{r} - \frac{\vec{r}\cdot\vec{r}}{r^3}\right) f' + \frac{\vec{r}\cdot\vec{r}}{r^2}f'' = \frac{n-1}{r}f' + f''$$

Therefore $\displaystyle \Delta f = \frac{1}{r^{n-1}}\frac{d}{dr}\left(r^{n-1}f'\right)$.


I don't know the history of writing the Laplacian of a radial scalar field this way. But it does make certain calculations fast. For example, in 3D do the scalar fields $f_1(r) = A/r$ or $f_2(r) = B/r^2$ have zero Laplacian?

It also makes calculating the Fourier transform of such a Laplacian easier.

Solution 2:

Here's an alternative, it uses some heavy machinery (if some points are unclear perhaps the comment at the end might help) but casts a little light on the symmetry of the situation.

Let's distinguish $f$ from $\phi$ with $f(x)=\phi(\|x\|)$, and from now on we fix $x$ and put $\|x\|=r$. $B_{r}$ is the ball centered in the origin with radius $r$

We have that \begin{equation*} s_{n}r^{n-1}\phi'(r)=\int_{\partial B_{r}} \nabla f . \vec u \end{equation*} where $\vec u$ is the normal vector to $\partial B_{r}$, and $s_{n}$ is the surface area of $\partial B_{1}$.

By divergence theorem \begin{equation*} \int_{\partial B_{r}} \nabla f . \vec u = \int_{B_{r}} \Delta f \end{equation*}

But \begin{equation*} \frac{\partial}{\partial r}\left(\int_{B_{r}} \Delta f\right) = \int_{\partial B_{r}}\Delta f \end{equation*}

So \begin{equation*} \frac{\partial}{\partial r}\left(s_{n}r^{n-1}\phi'(r)\right)=\int_{\partial B_{r}} \Delta f \end{equation*}

Since $\Delta f$ is also a radial function \begin{equation*} \frac{1}{s_{n}r^{n-1}} \int_{B_{r}} \Delta f= \Delta f(x) \end{equation*} which concludes our proof (the $s_{n}$ cancel out).

A first problem with this argument is that it makes use of the fact that $\nabla f(x) = \phi'(\|x\|)\frac{x}{\|x\|}$ and that $\nabla f$ is also a radial function. Proving this properly requires more or less as much calculations as computing directly the laplacian. But this properties can be easily seen when picturing a radial function. The gradient must be orthogonal to the plane tangent to $B_{r}$ in $x$ (because $f$ is constant on $B_{r}$), so it is colinear to $x$ and must have norm $\phi'(r)$. The laplacian can be expressed from the local mean of $f$, so it should by invariant by isometry by radiality of $f$ and hence only depend on $r$.

The advantage of this proof is that it gives a faint idea of why a $r^{n-1}$ goes in and out the expression of the laplacian: the radiality of $f$ makes it possible to compute the laplacian on the whole $\partial B_{r}$ at once and not only in one point, so the $r^{n-1}$ appears when taking means over the $\partial B_{r}$.