PLUTO  4.4-patch2
Computation of the Right Hand Side

The right hand side is constructed by adding, in succession, 1-D contributions coming from different sweeps:

\[ \vec{U}^{n+1}_{ijk} - \vec{U}^n_{\ivec} = \vec{\cal R}^{(x_1)}_{\ivec} + \vec{\cal R}^{(x_2)}_{\ivec} + \vec{\cal R}^{(x_3)}_{\ivec} \]

It is called once for every direction at every integration stage. The single contribution, say for the $ x_1 $ direction, is computed in a quasi-conservative way as

\[ \vec{\cal R}^{(x_1)}_{\ivec} = -\frac{\Delta t}{\dvol_i} \left( {\cal A}_{i+\HALF}\vec{\cal F}_{i+\HALF} - {\cal A}_{i-\HALF}\vec{\cal F}_{i-\HALF}\right) + \Delta t\left(\vec{\cal S}^{geo}_i + \vec{\cal S}^{body}_i + \vec{\cal S}^{8w}_i\right) +\vec{\cal S}^{eglm}_i \]

where the indices $ j $ and $ k $ stay fixed, $ \dvol_d $ is the cell volume, $ \vec{\cal F}^{\cal A} $ is the Riemann solver flux (not including pressure terms), $ {\cal A} $ is the interface area. The source terms $ \vec{\cal S}^{geo}_i, \vec{\cal S}^{body}_i, \vec{\cal S}^{8w}_d $ are, respectively, the geometrical, body force, Powell's and extended GLM source terms (only for the MHD and RMHD modules). Body forces are added by calling the AddBodyForce() function.

The pressure term is separately discretized as a gradient operator and therefore it does not contain area factors nor it appears in the geometrical source terms. For instance, while incrementing the component of $ \vec{\cal R} $ corresponding to the normal momentum $ m_{x_1} $ we add

\[ {\cal R}_{\{x_1\},\ivec}^{[m_{x_1}]} \quad -= \frac{\Delta t}{\Delta l_i} \left(p_{i+\HALF} - p_{i-\HALF}\right) \]

where $ \Delta l_i $ is a line element.

The computation of the right hand side closely reflects the nature of the divergence or gradient operators involved in the original equations. Indeed, for scalar quantities such as density or energy, one simply has

\[ \pd{q}{t} + \nabla\cdot\vec{F} = 0\quad\Longrightarrow\quad \left\{\begin{array}{c@{+}c@{+}c@{+}c@{=}cl} \DS \pd{q}{t} & \DS \pd{F_x}{x} & \DS \pd{F_y}{y} & \DS \pd{F_z}{z} & \DS\; 0 &\qquad\textrm{Cartesian} \\ \noalign{\medskip} \DS \pd{q}{t} & \DS \frac{1}{r}\pd{(rF_{r})}{r} & \DS \frac{1}{r}\pd{F_{\phi}}{\phi} & \DS \pd{F_z}{z} & \DS \;0 & \qquad\textrm{Polar} \\ \noalign{\medskip} \DS \pd{q}{t} & \DS \frac{1}{r^2}\pd{(r^2F_{r})}{r} & \DS \frac{1}{rs}\pd{(s F_{\theta})}{\theta} & \DS \frac{1}{rs}\pd{F_\phi}{\phi} & \DS\; 0 &\qquad\textrm{Spherical} \\ \noalign{\medskip} \end{array}\right. \]

For vector quantities such as momentum and magnetic field, we exploit the symmetric or antisymmetric properties of the corresponding flux tensor.

\[ \pd{\vec{B}}{t} + \nabla\cdot\vec{\Omega} = 0 \; \Longrightarrow \; \left\{\begin{array}{c@{+}c@{+}c@{+}c@{=}c} \DS \pd{B_{r}}{t} & \DS & \DS \frac{1}{rs} \pd{(s \Omega_{\theta r})}{\theta} & \DS \frac{1}{rs} \pd{\Omega_{\phi r}}{\phi} & \DS 0 \\ \noalign{\medskip} \DS \pd{B_{\theta}}{t} & \DS \frac{1}{r} \pd{(r\Omega_{r\theta})}{r} & \DS & \DS \frac{1}{rs}\pd{\Omega_{\phi\theta}}{\phi} & \DS 0 \\ \noalign{\medskip} \DS \pd{B_{\phi}}{t} & \DS \frac{1}{r} \pd{(r\Omega_{r\phi})}{r} & \DS \frac{1}{r} \pd{ \Omega_{\theta\phi}}{\theta} & \DS & \DS 0 \\ \noalign{\medskip} \end{array}\right. \]

The geometrical source terms (right hand sides in the previous equations) are evaluated by averaging interface values at the cell-centered.