I'm wondering if anyone can help me understand energy conservation when using conservation law methods (i.e. Riemann solver, High-Resolution Wave-Propagation Methods) with the addition of source terms:
$ \frac{\partial \vec{Q}}{\partial t} + \nabla \cdot \vec{F} = \vec{\psi} $,
where, $\vec{Q}$, and $\vec{F}$ are the conserved quantities and their corresponding fluxes. $\vec{\psi}$ are the non-hyperbolic source terms.
A 1D example of my Euler equation setup with an external gravitational field is:
$ \frac{\partial}{\partial t}\left[ \rho \right] + \frac{\partial}{\partial z}\left[ \rho w \right] = 0$
$ \frac{\partial}{\partial t}\left[ \rho w \right] + \frac{\partial}{\partial z}\left[ \rho w^2 + p \right] = -\rho g$
$\frac{\partial}{\partial t}\left[ e \right] + \frac{\partial}{\partial z}\left[ (e+p)w\right] = -\rho g w$
where, $t$ is time, $z$ is the spatial coordinate, $\rho$ is density, $w$ is fluid velocity, $e$ is the energy (internal + kinetic), $g$ is the acceleration due to gravity in the $z$ direction, and $p$ is the pressure. The system is closed by defining an equation of state:
$p = (\gamma -1)e_{i} $,
where $e_{i} = e - \frac{1}{2}\rho w^2$ is the internal energy, and $\gamma$ is the ratio of specific heats.
The system of equations is solved using a wave-propagation, Riemann solver with an F-Wave approach using the PyClaw software (http://www.clawpack.org/doc/pyclaw/). This means the non-hyperbolic gravity source terms are included in the flux differencing (rather than solved using some fractional step source splitting technique). The advantage of including the gravity source terms in the flux differencing is that in an equilibrium case, the source terms will exactly cancel the pressure gradient and maintain numerical equilibrium. This has proven to work very well.
My question is regarding the conservation of energy, particularly the total energy that also includes gravitational energy, $ e_{total} = e_i + e_{KE} + e_{gravitational}$. Is there any expectation that the total energy $e_{total}$ should be conserved as precisely as the "conserved" quantity energy $e$ when the hyperbolic terms are modified with source terms? In a test problem where an equilibrium case is perturbed, I see excellent conservation of $e$, but not great conservation of the total energy $e_{total}$. I'm not sure what to make of it.
I can provide a more detailed explanation of the equilibrium setup, and the perturbed case if that is helpful.
Thanks, Wes