### Is there a way to compute the "correct" distributional answer, -4*pi*dirac_delta(r), as the Laplacian of 1/r (in spherical coordinates on Euclidean space)?

**Background:** I'm putting together an EM notebook for educational purposes, and would like to be able to show demonstrations of Gauss's Law (for example). However, for completely understandable and obvious reasons, integrating `(1/r^2*r_hat).div()`

over the unit ball just gives `0`

. Computing the flux of `(1/r^2*r_hat)`

through the unit sphere gives the expected result, of course, but I'd hoped to directly demonstrate that it's proportional to the volume integral.