The Conservative Form and the Non-conservative Form

The difference between a conservative form and a non-conservative form lies in the splitting of derivatives. Taking the processing of the term $\rho\mathbf{g}$ as an example.
As is known, gravity $\mathbf{g}$ is the negative gradient of gravitational potential, i.e.

$$ \mathbf{g}=-\nabla{U}, \tag{1} $$

where $U=-\frac{GM}{r}$ is the gravitational potential.
In the Cartesian coordinate system, the equation can be rewritten as

$$ \mathbf{g}=-\frac{\partial{U}}{\partial{x}}\mathbf{i}-\frac{\partial{U}}{\partial{y}}\mathbf{j}-\frac{\partial{U}}{\partial{z}}\mathbf{k}, \tag{2} $$

which can be discretized by Finite Difference Method.
When multiplied by $\rho$, the equation becomes

$$ \rho\mathbf{g}=-\rho\frac{\partial{U}}{\partial{x}}\mathbf{i}-\rho\frac{\partial{U}}{\partial{y}}\mathbf{j}-\rho\frac{\partial{U}}{\partial{z}}\mathbf{k}, \tag{3} $$

where $\rho$ cannot be put into the derivative, since

$$ \frac{\partial{(\rho{}U)}}{\partial{x}} =\rho\frac{\partial{U}}{\partial{x}} +U\frac{\partial{\rho}}{\partial{x}}\ne{}\rho\frac{\partial{U}}{\partial{x}}. \tag{4} $$

However, in a conservative form, for the convenience of discretization, the variable $\Omega$ is introduced, and the equation is written as

$$ \rho\mathbf{g}=\nabla{(\rho\Omega)}=\frac{\partial{(\rho\Omega)}}{\partial{x}}\mathbf{i}+\frac{\partial{(\rho\Omega)}}{\partial{y}}\mathbf{j}+\frac{\partial{(\rho\Omega)}}{\partial{z}}\mathbf{k}, \tag{5} $$

which can form a telescoping series during the discretization, which means when the terms are added up over a grid, only the boundary terms would remain. For example

$$ \frac{(\rho\Omega)_1-(\rho\Omega)_2}{\Delta{}x}+\frac{(\rho\Omega)_2-(\rho\Omega)_3}{\Delta{}x}+\frac{(\rho\Omega)_3-(\rho\Omega)_4}{\Delta{}x}=\frac{(\rho\Omega)_1-(\rho\Omega)_4}{\Delta{}x} \tag{6} $$

This form is suitable for Finite Volume Method.

The Conservative Form of the Momentum Equations

According to the Reynolds transport theorem

$$ \frac{\mathrm{d}\mathbf{B}_{\mathrm{sys}}}{\mathrm{d}t}=\frac{\mathrm{d}}{\mathrm{d}t}\big(\iiint_{\mathrm{cv}}\boldsymbol{\beta}\rho\mathrm{d}V\big)+\iint_{\mathrm{cs}}\boldsymbol{\beta}\rho(\boldsymbol{\beta}\cdot\boldsymbol{n})\mathrm{d}A, \tag{7} $$

when $\mathbf{B}=m\mathbf{v}$ and $\boldsymbol{\beta}=\frac{\mathrm{d}\mathbf{B}}{\mathrm{d}m}=\mathbf{v}$, the equation becomes

$$ \frac{\mathrm{d}(m\mathbf{\mathbf{v}})_{\mathrm{sys}}}{\mathrm{d}t}=\frac{\mathrm{d}}{\mathrm{d}t}\big(\iiint_{\mathrm{cv}}\mathbf{v}\rho\mathrm{d}V\big)+\iint_{\mathrm{cs}}\mathbf{v}\rho(\mathbf{v}\cdot\boldsymbol{n})\mathrm{d}A. \tag{8} $$

According to the divergence theorem,

$$ \iint_{\mathrm{S}}\mathbf{n}\cdot\mathbf{U}\mathrm{d}S=\iiint_{\mathrm{V}}\nabla\cdot\mathbf{U}\mathrm{d}V, \tag{9} $$

equation $(8)$ becomes

$$ \begin{aligned} \frac{\mathrm{d}(m\mathbf{\mathbf{v}})_{\mathrm{sys}}}{\mathrm{d}t}&=\frac{\partial}{\partial{t}}\big(\iiint_{\mathrm{cv}}\mathbf{v}\rho\mathrm{d}V\big)+\iint_{\mathrm{cs}}\mathbf{v}\rho(\mathbf{v}\cdot\boldsymbol{n})\mathrm{d}A. \\ &=\frac{\partial}{\partial{t}}\big(\iiint_{\mathrm{cv}}\mathbf{v}\rho\mathrm{d}V\big)+\iiint_{\mathrm{cv}}\nabla\cdot(\rho\mathbf{vv})\mathrm{d}V \\ &=\iiint_{\mathrm{cv}}\big[\frac{\partial{(\rho\mathbf{v}})}{\partial{t}}+\nabla\cdot(\rho\mathbf{vv})\big]\mathrm{d}V \end{aligned} \tag{10} $$

To keep the conservative form of the derivatives, one should not expand the term $\nabla\cdot(\rho\mathbf{v}\mathbf{v})$.

So, the conservative form of the momentum equation is written as

$$ \left\{ \begin{split} \quad\frac{\partial{(\rho\mathbf{v})}}{\partial{t}}+\nabla\cdot(\rho{}\mathbf{v}\mathbf{v})&=-\nabla{p}+\rho\mathbf{g}+\nabla{}\cdot\mathbf{T} \\ \rho\mathbf{g}&=\frac{\partial{(\rho\Omega)}}{\partial{x}}\mathbf{i}+\frac{\partial{(\rho\Omega)}}{\partial{y}}\mathbf{j}+\frac{\partial{(\rho\Omega)}}{\partial{z}}\mathbf{k}, \end{split} \right. \tag{11} $$

The Dyadic Product of $\mathbf{vv}$

The operation of $\mathbf{vv}$ is neither a cross product or a dot product. It is called a dyadic product, which can also be written as $\mathbf{v\otimes{}v}$.

$$ \left. \begin{split} \mathbf{vv} &= (v_x\mathbf{i} + v_y\mathbf{j} + v_z\mathbf{k} )(v_x\mathbf{i} + v_y\mathbf{j} + v_z\mathbf{k} )\\ &=\mathbf{ii} v_xv_x+ \mathbf{ij} v_xv_y + \mathbf{ik} v_xv_z \\ &\quad+\mathbf{ji} v_yv_x+ \mathbf{jj} v_yv_y + \mathbf{jk} v_yv_z \\ &\quad+\mathbf{ki} v_zv_x+ \mathbf{kj} v_zv_y + \mathbf{kk} v_zv_z \end{split} \right\} \Rightarrow \mathbf{vv} = \begin{pmatrix} v_xv_x & v_xv_y & v_xv_z \\ v_yv_x & v_yv_y & v_yv_z \\ v_zv_x & v_zv_y & v_zv_z \\ \end{pmatrix} \tag{12} $$

$$ \left. \begin{split} \nabla \cdot{}(\mathbf{\rho{}vv}) &= (\frac{\partial}{\partial{x}}\mathbf{i}+\frac{\partial}{\partial{y}}\mathbf{j}+\frac{\partial}{\partial{z}}\mathbf{k})\cdot{}(v_x\mathbf{i} + v_y\mathbf{j} + v_z\mathbf{k} )(v_x\mathbf{i} + v_y\mathbf{j} + v_z\mathbf{k} )\\ &=(\frac{\partial}{\partial{x}}\mathbf{i}+\frac{\partial}{\partial{y}}\mathbf{j}+\frac{\partial}{\partial{z}}\mathbf{k})\cdot(\mathbf{ii} v_xv_x+ \mathbf{ij} v_xv_y + \mathbf{ik} v_xv_z \\ &\quad+\mathbf{ji} v_yv_x+ \mathbf{jj} v_yv_y + \mathbf{jk} v_yv_z \\ &\quad+\mathbf{ki} v_zv_x+ \mathbf{kj} v_zv_y + \mathbf{kk} v_zv_z) \\ &= \frac{\partial{v_xv_x}}{\partial{x}}\mathbf{i} + \frac{\partial{v_xv_y}}{\partial{x}}\mathbf{j} + \frac{\partial{v_xv_z}}{\partial{x}}\mathbf{k} \\ &\quad+ \frac{\partial{v_yv_x}}{\partial{y}}\mathbf{i} + \frac{\partial{v_yv_y}}{\partial{y}}\mathbf{j} + \frac{\partial{v_yv_z}}{\partial{y}}\mathbf{k} \\ &\quad+ \frac{\partial{v_zv_x}}{\partial{x}}\mathbf{i} + \frac{\partial{v_zv_y}}{\partial{x}}\mathbf{j} + \frac{\partial{v_zv_z}}{\partial{x}}\mathbf{k} \\ \end{split} \right\} \Rightarrow \nabla \cdot{}(\mathbf{\rho{}vv}) = \begin{pmatrix} \frac{\partial{v_xv_x}}{\partial{x}} + \frac{\partial{v_yv_x}}{\partial{y}} + \frac{\partial{v_zv_x}}{\partial{z}} \\ \frac{\partial{v_xv_y}}{\partial{x}} + \frac{\partial{v_yv_y}}{\partial{y}} + \frac{\partial{v_zv_y}}{\partial{z}} \\ \frac{\partial{v_xv_z}}{\partial{x}} + \frac{\partial{v_yv_z}}{\partial{y}} + \frac{\partial{v_zv_z}}{\partial{z}} \end{pmatrix} \tag{13} $$

Where, for example, $\mathbf{i\cdot{}ij } = \mathbf{(i\cdot{}i)j } = \mathbf{j } $ and $\mathbf{i\cdot{}ji } = \mathbf{(i\cdot{}j)i } = \mathbf{0} $.

Simmilarly,

$$ \left. \begin{split} \nabla \cdot{}\mathbf{T} &=(\frac{\partial}{\partial{x}}\mathbf{i}+\frac{\partial}{\partial{y}}\mathbf{j}+\frac{\partial}{\partial{z}}\mathbf{k})\cdot(\mathbf{ii} \tau_{xx}+ \mathbf{ij} \tau_{xy} + \mathbf{ik} \tau_{xz} \\ &\quad+\mathbf{ji} \tau_{yx}+ \mathbf{jj} \tau_{yy} + \mathbf{jk} \tau_{yz} \\ &\quad+\mathbf{ki} \tau_{zx}+ \mathbf{kj} \tau_{zx} + \mathbf{kk} \tau_{zz}) \\ &= \frac{\partial{\tau_{xx}}}{\partial{x}}\mathbf{i} + \frac{\partial{\tau_{xy}}}{\partial{x}}\mathbf{j} + \frac{\partial{\tau_{xz}}}{\partial{x}}\mathbf{k} \\ &\quad+ \frac{\partial{\tau_{yx}}}{\partial{y}}\mathbf{i} + \frac{\partial{\tau_{yy}}}{\partial{y}}\mathbf{j} + \frac{\partial{\tau_{yz}}}{\partial{y}}\mathbf{k} \\ &\quad+ \frac{\partial{\tau_{zx}}}{\partial{x}}\mathbf{i} + \frac{\partial{\tau_{zy}}}{\partial{x}}\mathbf{j} + \frac{\partial{\tau_{zz}}}{\partial{x}}\mathbf{k} \\ \end{split} \right\} \Rightarrow \nabla \cdot{} \mathbf{T}= \begin{pmatrix} \frac{\partial{\tau_{xx}}}{\partial{x}} + \frac{\partial{\tau_{yx}}}{\partial{y}} + \frac{\partial{\tau_{zx}}}{\partial{z}} \\ \frac{\partial{\tau_{xy}}}{\partial{x}} + \frac{\partial{\tau_{yy}}}{\partial{y}} + \frac{\partial{\tau_{zy}}}{\partial{z}} \\ \frac{\partial{\tau_{xz}}}{\partial{x}} + \frac{\partial{\tau_{yz}}}{\partial{y}} + \frac{\partial{\tau_{zz}}}{\partial{z}} \end{pmatrix} \tag{14} $$

The Matrix Form of the Momentum Equations

Expand the momentum equation $(11)$
and get

$$ \begin{aligned} &\frac{\partial{}}{\partial{t}} \begin{pmatrix}\rho{}v_x\\\rho{}v_y\\\rho{}v_z\end{pmatrix} +\frac{\partial{}}{\partial{x}} \begin{pmatrix}\rho{}v_x^2\\\rho{}v_xv_y\\\rho{}v_xv_z\end{pmatrix} +\frac{\partial{}}{\partial{y}} \begin{pmatrix}\rho{}v_yv_x\\\rho{}v_y^2\\\rho{}v_yv_z\end{pmatrix} +\frac{\partial{}}{\partial{z}} \begin{pmatrix}\rho{}v_zv_x\\\rho{}v_zv_y\\\rho{}v_z^2\end{pmatrix} \\ =& -\begin{pmatrix}\frac{\partial{p}}{\partial{x}}\\\frac{\partial{p}}{\partial{y}}\\\frac{\partial{p}}{\partial{z}}\end{pmatrix} + \begin{pmatrix}\frac{\partial{(\rho\Omega)}}{\partial{x}}\\\frac{\partial{(\rho\Omega)}}{\partial{y}}\\\frac{\partial{(\rho\Omega)}}{\partial{z}}\end{pmatrix} + \begin{pmatrix}\frac{\partial{\tau_{xx}}}{\partial{x}}+\frac{\partial{\tau_{yx}}}{\partial{y}}+\frac{\partial{\tau_{zx}}}{\partial{z}} \\ \frac{\partial{\tau_{xy}}}{\partial{x}}+\frac{\partial{\tau_{yy}}}{\partial{y}}+\frac{\partial{\tau_{zy}}}{\partial{z}} \\ \frac{\partial{\tau_{xz}}}{\partial{x}}+\frac{\partial{\tau_{yz}}}{\partial{y}}+\frac{\partial{\tau_{zz}}}{\partial{z}}\end{pmatrix}, \end{aligned} \tag{15} $$

which can then be processed into

$$ \frac{\partial{}}{\partial{t}} \begin{pmatrix}\rho{}v_x\\\rho{}v_y\\\rho{}v_z\end{pmatrix} +\frac{\partial{}}{\partial{x}} \begin{pmatrix}\rho{}(v_x^2-\Omega)+p-\tau_{xx}\\\rho{}v_xv_y-\tau_{xy}\\\rho{}v_xv_z-\tau_{xz}\end{pmatrix} +\frac{\partial{}}{\partial{y}} \begin{pmatrix}\rho{}v_yv_x-\tau_{yx}\\\rho{}(v_y^2-\Omega)+p-\tau_{yy}\\\rho{}v_yv_z-\tau_{yz}\end{pmatrix} +\frac{\partial{}}{\partial{z}} \begin{pmatrix}\rho{}v_zv_x-\tau_{zx}\\\rho{}v_zv_y-\tau_{zy}\\\rho{}(v_z^2-\Omega)+p-\tau_{zz}\end{pmatrix}=0 \tag{16} $$