Solution (source code)

= Solution

Use the <BGK model> from <1/b>. Multiply it by $v_i$ and integrate over <velocity space>. Matching the <collision invariants> makes the collision contribution vanish. For the external <force> term, <integration by parts> gives
$$
-\partial_jU\int v_i\partial_{v_j}f\,d^3v=\rho\,\partial_iU.
$$
Define the <momentum flux> and <force density> by
$$
\Pi_{ij}=\int v_iv_jf\,d^3v,\qquad F_i=-\rho\,\partial_iU.
$$
The resulting <momentum conservation> equation is
$$
\boxed{\partial_t(\rho u_i)+\partial_j\Pi_{ij}=F_i,\qquad \boldsymbol\Pi=\rho\mathbf u\mathbf u+\mathbf P,\qquad \mathbf F=-\rho\nabla U.}
$$
The <pressure tensor> is the second central <velocity moment>,
$$
P_{ij}=\int(v_i-u_i)(v_j-u_j)f\,d^3v.
$$
The zeroth <velocity moment> also supplies the <continuity equation> $\partial_t\rho+\partial_j(\rho u_j)=0$, which will be needed in <1/d>.