Divergence theorem and singularities Consider a function $\overline{f}=\frac{1}{r^2}\widehat{r}$ where $\widehat{r}$ is the unit vector in the radial direction. Calculate the divergence of this function over a sphere of radius $R$ centered at the origin.
The divergence theorem states:
$\iiint_V \nabla \cdot  \overline{f}\,dV\,=\,\unicode{x222F}_{S(V)}\overline{f}\cdot\,\widehat{n}\,\,dS$
Solving the RHS,
$\unicode{x222F}_{S(V)}\overline{f}\cdot\,\widehat{n}\,\,dS = \unicode{x222F}_{S(V)} (\frac{1}{r^2} )\cdot(r^2.sin\theta .d\theta .d\phi)=4\pi$
Solving the LHS,
In the spherical co-ordinate system we have: $\nabla \cdot \overline{f}=\frac{1}{r^2}.\frac{\partial(r^2.f_{r})}{\partial r}$   (Considering only $\widehat{r}$)
But,
$\frac{\partial(r^2.\frac{1}{r^2})}{\partial r}=\frac{\partial (1)}{\partial{r}}=0$
$\therefore \,\,\,\,\iiint_V \nabla \cdot  \overline{f}\,dV=0$
I understand this has something to do with the singularity at the origin. But I don't know what to do in this case. Is $4\pi$ a valid answer? Can someone explain the correct way to solve this discrepancy?
 A: The divergence of $ \overline{f} $ is zero everywhere except at the origin where it is undefined, because of the $ \frac{1}{r^2} $ term.
However, the integral of $ \overline{f} $ over $V$ is not zero, because the origin acts as a "source of divergence". This is analogous to the Dirac delta function, which is zero everywhere but at the origin, and whose integral is 1.
Note that since the origin is the only "source of divergence", the flux of $ \overline{f} $ through the sphere has the same value for every radius $R$ — the value of $ 4\pi $ is correct.
If the domain of integration didn't include the origin, this value would have been zero.
A: You can apply the divergence theorem only to continuously differentiable vector fields, so (as you hint) you can't apply it in the unit ball volume. What you can do is apply it on any $V = B(1) \setminus B(a)$, where $B(a)$ is the ball of radius $a$ centred on the origin, and take the limit of $a$ going to zero.
Then you'll find that the outward fluxes on $S(a)=\partial B(a)$ and $S(1)$ are equal and opposite, $-4\pi$ and $4\pi$.
This, by the way, is also the technique to calculate the distribution that corresponds to the divergence of your vector field, which is a Dirac mass at the origin.
