Solution (source code)

= Solution

Multiply the wave equation by $\rho_0\phi_t/c_0^2$ and use product rules. This gives
$$
\frac{\partial}{\partial t}
\left[\frac{\rho_0}{2}|\nabla\phi|^2
+\frac{\rho_0}{2c_0^2}\phi_t^2\right]
+\nabla\cdot(-\rho_0\phi_t\nabla\phi)=0.
$$
Therefore
$$
\boxed{
\mathcal E=\frac{\rho_0}{2}|\nabla\phi|^2
+\frac{\rho_0}{2c_0^2}\phi_t^2
=\frac{\rho_0}{2}|\mathbf u|^2+\frac{p'^2}{2\rho_0c_0^2},
\qquad
\mathbf I=-\rho_0\phi_t\nabla\phi=p'\mathbf u.
}
$$
These are respectively the kinetic-plus-compressional acoustic energy density and its flux.