Let $U \subset \mathbb{R}^3$ be a bounded subset with smooth boundary.
Let $Y \colon \mathbb{R}^3 \to \mathbb{R}^3 $ be a smooth vector field. I know that the first variation of the volume of $U$ w.r.t. $Y$ is given by $$ \delta_Y|U| = \int_U \text{div}Y \, dx_1 dx_2dx_3. $$ I know that it follows from some standard computation, but I can't remember how to derive it.
Let $\{\phi_t\}_{t\in (-\epsilon, \epsilon)}$ be the flow corresponding to $Y$. Then
\begin{align} \delta_Y [U] &= \frac{d}{dt}\bigg|_{t=0} \int_{\phi_t(U)} \mathrm d x_1\mathrm d x_2\mathrm d x_3\\ &= \frac{d}{dt}\bigg|_{t=0} \int_{U} \phi_t^*(\mathrm d x_1\mathrm d x_2\mathrm d x_3) \\ &= \frac{d}{dt}\bigg|_{t=0} \int_U \mathrm d\phi_1 \mathrm d \phi_2 \mathrm d \phi_3 \\ &= \int_U \mathrm dY_1 \mathrm d x_2 \mathrm dx_3+\mathrm dx_1 \mathrm d Y_2 \mathrm dx_3+\mathrm dx_1 \mathrm d x_2 \mathrm dY_3 \\ &= \int_U \left(\frac{\partial Y_1}{\partial x_1}+\frac{\partial Y_2}{\partial x_2}+\frac{\partial Y_3}{\partial x_3}\right) \mathrm dx_1\mathrm dx_2\mathrm dx_3\\ &=\int_U \text{div} Y\mathrm dx_1\mathrm dx_2\mathrm dx_3. \end{align}
Note we write $\phi_t =(\phi_1, \phi_2, \phi_3)$ and $Y=(Y_1, Y_2, Y_3)$.