I was trying to derive $\nabla \times \vec{B}\left(\vec{r}\right) = \mu_0 \vec{J}\left(\vec{r}\right)$ from: $$ \vec{B} \left( \vec{r} \right) = \frac{\mu_0}{4\pi} \iiint_{V_s} \frac{\vec{J}\left(\vec{r}_s\right)}{\left|\vec{r}-\vec{r}_s\right|^3} \times \left(\vec{r}-\vec{r}_s\right) \space dV\left(\vec{r}_s\right) $$ So I first substitute $\nabla\frac{1}{|\vec{r}-\vec{a}|} = -\frac{\vec{r}-\vec{a}}{|\vec{r}-\vec{a}|^3}$ and get: $$ \vec{B} \left( \vec{r} \right) = -\frac{\mu_0}{4\pi} \iiint_{V_s} \vec{J}\left(\vec{r}_s\right) \times \nabla_{\vec{r}} \frac{1}{\left|\vec{r}-\vec{r}_s\right|} \space dV\left(\vec{r}_s\right) $$ Which, because $\vec{J}\left(\vec{r}_s\right)$ doesn't depend on $\vec{r}$, is equal to: $$ \vec{B} \left( \vec{r} \right) = \frac{\mu_0}{4\pi} \iiint_{V_s} \nabla_{\vec{r}} \times \frac{\vec{J}\left(\vec{r}_s\right)}{\left|\vec{r}-\vec{r}_s\right|} \space dV\left(\vec{r}_s\right) $$ Now if I take the curl of both side, I get: $$ \nabla_{\vec{r}} \times \vec{B} \left( \vec{r} \right) = \frac{\mu_0}{4\pi} \iiint_{V_s} \nabla_{\vec{r}} \times \left( \nabla_{\vec{r}} \times \frac{\vec{J}\left(\vec{r}_s\right)}{\left|\vec{r}-\vec{r}_s\right|} \right) \space dV\left(\vec{r}_s\right) $$ Which is: $$ \frac{\mu_0}{4\pi} \iiint_{V_s} \left( \nabla_{\vec{r}}\left(\nabla_{\vec{r}} \cdot \frac{\vec{J}\left(\vec{r}_s\right)}{\left|\vec{r}-\vec{r}_s\right|} \right) - {\nabla^2}_{\vec{r}} \frac{\vec{J}\left(\vec{r}_s\right)}{\left|\vec{r}-\vec{r}_s\right|} \right) \space dV\left(\vec{r}_s\right) $$ Which, because $\vec{J}\left(\vec{r}_s\right)$ doesn't depend on $\vec{r}$, is: $$ \frac{\mu_0}{4\pi} \iiint_{V_s} \left( \nabla_{\vec{r}} \left( \vec{J}\left(\vec{r}_s\right) \cdot \nabla_{\vec{r}} \frac {1} {\left|\vec{r}-\vec{r}_s\right|} \right) + \vec{J}\left(\vec{r}_s\right) \space 4\pi \delta^3 \left(\vec{r}-\vec{r}_s\right) \right) \space dV\left(\vec{r}_s\right) $$ Which is: $$ \left[ \frac{\mu_0}{4\pi} \iiint_{V_s} \nabla_{\vec{r}} \left( \vec{J}\left(\vec{r}_s\right) \cdot \nabla_{\vec{r}} \frac {1} {\left|\vec{r}-\vec{r}_s\right|} \right) \space dV\left(\vec{r}_s\right) \right] \Large{ + \space\mu_0 \vec{J}\left(\vec{r}\right)} $$ Which is: $$ \frac{\mu_0}{4\pi} \left[ \nabla_{\vec{r}} \iiint_{V_s} \left( \vec{J}\left(\vec{r}_s\right) \cdot \nabla_{\vec{r}} \frac {1} {\left|\vec{r}-\vec{r}_s\right|} \right) \space dV\left(\vec{r}_s\right) \right] \Large{+ \space\mu_0 \vec{J}\left(\vec{r}\right)} $$ Now my question is, given $\nabla \cdot \vec{J}\left(\vec{r}_s\right) = 0$, how can I show that $ \nabla \iiint_{V_s} \left( \vec{J}\left(\vec{r}_s\right) \cdot \nabla \frac {1} {\left|\vec{r}-\vec{r}_s\right|} \right) \space dV\left(\vec{r}_s\right) $ reduces to zero? Which expression reduces to $\nabla \cdot \vec{J}\left(\vec{r}_s\right)$, and how? What is the applicable identity?
Attempt
So to attack the integral $ \iiint_{V_s} \left( \vec{J}\left(\vec{r}_s\right) \cdot \nabla \frac {1} {\left|\vec{r}-\vec{r}_s\right|} \right) \space dV\left(\vec{r}_s\right) $, I first switch the variable over which the gradient is taken from $\vec{r}$ to $\vec{r}_s$: $$ \nabla_{\vec{r}_s}\frac {1} {\left|\vec{r}-\vec{r}_s\right|} = -\nabla_{\vec{r}}\frac {1} {\left|\vec{r}-\vec{r}_s\right|} $$ My integral, therefore, becomes: $$ -\iiint_{V_s} \left( \vec{J}\left(\vec{r}_s\right) \cdot \nabla_{\vec{r}_s} \frac {1} {\left|\vec{r}-\vec{r}_s\right|} \right) \space dV\left(\vec{r}_s\right) $$ Which is, doing a vector integration by parts: $$ \iiint_{V_s} \frac {\nabla_{\vec{r}_s} \cdot \vec{J}\left(\vec{r}_s\right)} {\left|\vec{r}-\vec{r}_s\right|} dV\left(\vec{r}_s\right) - \oint_{\partial V_s} \frac {\vec{J}\left(\vec{r}_s\right)} {\left|\vec{r}-\vec{r}_s\right|} \cdot d\vec{S}\left(\vec{r}_s\right) $$ The first term is zero because we know that $\nabla\cdot\vec{J}=0$.
I am guessing that the second term also zero because, by definition, there is no current density outside of $V_s$. (This is correct, as confirmed by the answer below)
You're basically right. More precisely, if you're integrating over a finite volume, you've already implicitly assumed that the current density outside of it vanishes; otherwise it would be wrong to calculate the magnetic field by integrating only over that volume. If you're integrating over all space, there's no "outside", and the corresponding assumption, which is usually made in electrodynamics, is that the the sources and fields decay sufficiently quickly towards infinity. In this case, the current has to decay faster than $1/r$ for the boundary term to vanish in the limit.
Note that you can find this derivation in Wikipedia and on p. $137$ of archive.org's version of Jackson's Classical Electrodynamics (which is nice, since $137$ is almost exactly the reciprocal of the fine-structure constant :-).