A comprehensive discussion of SPH can be found in the reviews of Hernquist and Katz , Benz , and Monaghan [198, 199]. The non-relativistic SPH equations are briefly discussed in Section 9.6. The capabilities and limits of SPH are explored, e.g., in [269, 16, 167, 275], and the stability of the SPH algorithm is investigated in .
The SPH equations for special relativistic flows have been first formulated by Monaghan . Monaghan and Price  showed how the equations of motion for the SPH method may be derived from a variational principle for both non-relativistic and (special and general) relativistic flows when there is no dissipation. For relativistic flows the SPH equations given in Section 9.6 can be used except that each SPH particle carries baryons instead of mass [198, 53]. Hence, the rest mass of particle is given by , where is the baryon rest mass (if the fluid is made of baryons). Transforming the notation used in  to ours, the continuity equation, the momentum, and the total energy equations for particle are given by (unit of velocity is c)
Special relativistic flow problems have been simulated with SPH by [151, 134, 172, 174, 53, 262]. Extensions of SPH capable of treating general relativistic flows have been considered by [134, 150, 262, 202, 204]. Concerning relativistic SPH codes the artificial viscosity is the most critical issue. It is required to handle shock waves properly, and ideally it should be predicted by a relativistic kinetic theory for the fluid. However, unlike its Newtonian analogue, the relativistic theory has not yet been developed to the degree required to achieve this.
For Newtonian SPH, Lattanzio et al.  have shown that a viscosity quadratic in the velocity divergence is necessary in high Mach number flows. They proposed a form such that the viscous pressure could be simply added to the fluid pressure in the equation of motion and the energy equation. As this simple form of the artificial viscosity has known limitations, they also proposed a more sophisticated form of the artificial viscosity terms, which leads to a modified equation of motion. This artificial viscosity works much better, but it cannot be generalized to the relativistic case in a consistent way. Utilizing an equation for the specific internal energy, both Mann  and Laguna et al.  use such an inconsistent formulation. Their artificial viscosity term is not included in the expression of the specific relativistic enthalpy. In a second approach, Mann  allows for a time-dependent smoothing length and SPH particle mass, and further proposes an SPH variant based on the total energy equation. Lahy  and Siegler and Riffert  use a consistent artificial viscosity pressure added to the fluid pressure. Siegler and Riffert  have also formulated the hydrodynamic equations in conservation form (see also ).
Monaghan  incorporates concepts from Riemann solvers into SPH (see also ). For this reason he also proposes to use a total energy equation in SPH simulation instead of the commonly used internal energy equation, which would involve time derivatives of the Lorentz factor in the relativistic case. Chow and Monaghan  have extended this concept and have proposed an SPH algorithm, which gives good results when simulating an ultra-relativistic gas. In both cases the intention was not to introduce Riemann solvers into the SPH algorithm, but to use them as a guide to improve the artificial viscosity required in SPH. Multi-dimensional simulations of general relativistic flows (in a given time-independent metric) using the SPH formulation of Monaghan and Price  and the SPH algorithm of Chow and Monaghan  have been performed by Muir .
In Roe’s Riemann solver , as well as in its relativistic variant proposed by Eulerdink [83, 84] (see Section 3.4), the numerical flux is computed by solving a locally linear system, and depends on both the eigenvalues and (left and right) eigenvectors of the Jacobian matrix associated to the fluxes and on the jumps in the conserved physical variables (see Equations (36) and (37)). Monaghan  realized that an appropriate form of the dissipative terms and for the interaction between particles a and b can be obtained by treating the particles as the equivalent of left and right states taken with reference to the line joining the particles. The quantity corresponding to the eigenvalues (wave propagation speeds) is an appropriate signal velocity (see below), and that equivalent to the jump across characteristics is a jump in the relevant physical variable. For the artificial viscosity tensor, , Monaghan  assumes that the jump in velocity across characteristics can be replaced by the velocity difference between a and b along the line joining them.
With these considerations in mind, Chow and Monaghan  proposed for in the relativistic case the forma and b are approaching, and = 0 otherwise. Here K = 0.5 is a dimensionless parameter, which is chosen to have the same value as in the non-relativistic case . is the average baryon number density, which has to be present in Equation (62), because the pressure terms in the summation of Equation (101) (see Section 9.6) have an extra density in the denominator arising from the SPH interpolation. Furthermore, b to a, and .
The dissipation term in the energy equation is derived in a similar way, and is given by a and b are approaching, and otherwise. involves the energy , which is identical to (see Equation (61)) except that is replaced by .
To determine the signal velocity, Chow and Monaghan  (and Monaghan  in the non-relativistic case) start from the (local) eigenvalues, and hence the wave velocities and of one-dimensional relativistic hydrodynamic flows. Again considering particles a and b as the left and right states of a Riemann problem with respect to motions along the line joining the particles, the appropriate signal velocity is the speed of approach (as seen in the computing frame) of the signal sent from a towards b and that from b to a. This is the natural speed for the sharing of physical quantities, because when information about the two states meets it is time to construct a new state. This speed of approach should be used when determining the size of the time step by the Courant condition (for further details see ).
Chow and Monaghan  have demonstrated the performance of their Riemann problem guided relativistic SPH algorithm by calculating several shock tube problems involving ultra-relativistic speeds up to v = 0.9999. The algorithm gives good results, but finite volume schemes based on Riemann solvers give more accurate results and can handle even larger speeds (see Section 6).
© Max Planck Society and the author(s)