(-\beta^k S_i + \alpha  S^k_{\ i}  )] =
    \sqrt{\gamma} \left[\frac{\alpha}{2}
  S^{lm} \partial_i \gamma_{lm}+ S_k \partial_i \beta^k 
- (\tau +
  D) \partial_i \alpha \right]\,. 
\label{momentum_evol}
\end{split}\end{aligned}$$

\[eq:rmhdconservative\]

The fluid variables $D,~U$, $S_{i}$, and $S_{ij}$ are the conserved rest-mass density, the conserved energy density, the conserved momentum, and the fully spatial projection of the energy-momentum tensor, respectively. Their explicit definitions are therefore $$\begin{aligned}
\label{Tmunu_decomposition2}
   D &:=& \rho W \,, \\
   \tau &:=& U - D = \rho h W^2 - p + \frac{1}{2} (E^2 + B^2) - D\,, \\
   S_{i} &:=& \rho h W^2 v_{i} + \epsilon_{ijk} E^j B^k\,, \\
   S_{ij} &:=& \rho h W^2 v_{i} v_{j} + \gamma_{ij} p 
 -E_i E_j - B_i B_j + \frac{1}{2} \gamma_{ij} (E^2 + B^2)\,.\end{aligned}$$

Here, $W=\alpha u^0 = u_i/v_i$ is the Lorentz factor, where $u^{\mu}$ are the components of the fluid 4-velocity and $v^i$ are the components of the 3-velocity as measured by the Eulerian observer. Furthermore, $h:=(e+p)/\rho = 1+\epsilon + p/\rho$ is the enthalpy, with $e = \rho
(1+\epsilon)$ the total energy density, $p$ the pressure, $\epsilon$ the specific internal energy, and $\rho$ the rest-mass density [@Rezzolla_book:2013]. The important difference between the RMHD and IMHD equations is that they involve stiff relaxation terms that 