In my vector analysis class, we were given a vector field
$\mathbf{F}(\mathbf{r})= \frac{1}{r^3}\mathbf{r}$
This vector field's divergence is 0, so according to the divergence theorem the flux over any closed object's surface $S_1$ within this field is 0. However, when the object has two surfaces, one outer surface $S_1$ and an inner surface, shaped like a ball with radius $\epsilon$, centered at $(0,0)$, when we calculate the flux over the surface $S_1$, it comes out as $4\pi$.
The exercise itself didn't require us to answer why or how, but I can't understand the physical interpretation of what's going on. Why does adding a second surface, $S_\epsilon$, to the object change the flux over the original surface, $S_1$?
Let $E$ be the region of space outside $S_{\epsilon}$ and inside $S_1$.
Since $F$ is continuous on $E$ which has boundary $S_1 \cup S_{\epsilon}$ we have by divergence theorem $$\iiint_{E}\text{div}(F)\mathrm{d}V=\iint_{S_1}\left(F\cdot n\right)\mathrm{d}S-\iint_{S_{\epsilon}}(F\cdot n)\mathrm{d}S$$ where $n$ is normal to the surface that points away from the origin. Because $\text{div}(F)=0$, $$\iint_{S_1}\left(F\cdot n\right)\mathrm{d}S=\iint_{S_{\epsilon}}(F\cdot n)\mathrm{d}S$$ Now $n=\frac{(x,y,z)}{\epsilon}$ for $\iint_{S_{\epsilon}}(F\cdot n)\mathrm{d}S$ and $x^2+y^2+z^2=\epsilon^2$ for any $(x,y,z)\in S_{\epsilon}$ so that $$\begin{eqnarray*}\iint_{S_\epsilon}(F\cdot n)\mathrm{d}S &=& \iint_{S_\epsilon}\left(\frac{(x,y,z)}{\epsilon^3} \cdot \frac{(x,y,z)}{\epsilon}\right)\mathrm{d}S \\ &=& \frac{1}{\epsilon^2}\iint_{S_{\epsilon}}\mathrm{d}S \\ &=& \frac{4\pi \epsilon^2}{\epsilon^2} \\ &=& 4\pi\end{eqnarray*}$$ We conclude $\iint_{S_1}(F\cdot n)\mathrm{d}S=4\pi$ as well.