In THIS ANSWER, I provided a rigorous development to show that the Dirac Delta can be regularized as
$$\bbox[5px,border:2px solid #C0A000]{\lim_{a\to 0} \nabla \cdot \vec \psi(\vec r;a)=4\pi \delta(\vec r)}$$
where
$$\vec \psi(\vec r;a)=\frac{\vec r}{(r^2+a^2)^{3/2}} $$
Here, we will use a formal approach that can easily be made rigorous as in the aforementioned post.
The interpretation of the regularization of the Dirac Delta given by
$$\delta(\vec r-\vec a)\sim \lim_{\alpha \to 0}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}$$
is that for any smooth test function $f(\vec r)$, we have
$$\lim_{\alpha \to 0}\int_{\mathscr{R}^3}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}\,f(\vec r)\,dV=f(\vec a)$$
We split the integral into the sum of two integrals
$$\begin{align}
\int_{\mathscr{R}^3}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}\,f(\vec r)\,dV&=\int_{\mathscr{R}^3-B_{\vec a}(\delta)}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}\,f(\vec r)\,dV\\\\
&+\int_{B_{\vec a}(\delta)}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}\,f(\vec r)\,dV
\end{align}$$
where $B_{\vec a}$ is a sphere of radius $\delta$ centered at $\vec a$. Note that as $\alpha \to 0$, the first integral vanishes. Heuristically, for small $\delta$, the second integral can be approximated by exploiting the continuity of $f$. Then, we can write
$$\begin{align}
\lim_{\alpha \to 0}\int_{\mathscr{R}^3}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}\,f(\vec r)\,dV&=f(\vec a)\lim_{\alpha \to 0}\int_{B_{\vec a}(\delta)}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}\,dV\\\\
&=f(\vec a)\lim_{\alpha \to 0}\frac{\alpha}{\pi^2}\int_0^{2\pi}\int_0^\pi \int_0^{\delta}\frac{1}{(\alpha^2+r^2)^2}r^2\sin(\theta)\,dr\,d\theta\,d\phi\\\\
&=f(\vec a)\lim_{\alpha \to 0}\left(\frac{\alpha}{\pi^2}\,(2\pi)\,(2)\,\frac{\pi}{4\alpha}\right)\\\\
&=f(\vec a)
\end{align}$$
as was to be shown!
Therefore, we have for any smooth test function $f(\vec r)$
$$\bbox[5px,border:2px solid #C0A000]{\lim_{\alpha \to 0}\int_{\mathscr{R}^3}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}\,f(\vec r)\,dV=f(\vec a)}$$
and it is in this sense that
$$\bbox[5px,border:2px solid #C0A000]{\delta(\vec r-\vec a)\sim \lim_{\alpha \to 0}\frac{\alpha /\pi^2}{(\alpha^2+|\vec r-\vec a|^2)^2}}$$