PLUTO
4.4-patch2
|
Implementation of the two-shock Riemann solver for the relativistic hydro equations. More...
#include "pluto.h"
Functions | |
static double | TwoShock_Lorentz (double *U, int n) |
static void | TwoShock_Shock (double, double, double, double, double, double, double, double *, double *, double *, int) |
void | TwoShock_Solver (const Sweep *sweep, int beg, int end, double *cmax, Grid *grid) |
double | TwoShock_RarefactionSpeed (double w[], int iside) |
Solve the Riemann problem for the relativistic hydrodynamics equations using the two-shock approach described by Mignone, Plewa and Bodo (2005). The formulation works with IDEAL or TAUB equation of state.
On input, this function takes left and right primitive state vectors stateL->v
and stateR->v
at zone edge i+1/2
. On output, return flux and pressure vectors at the same interface i+1/2
(note that the i
refers to i+1/2
).
Also during this step, compute maximum wave propagation speed (cmax) for explicit time step computation.
Reference:
|
static |
Compute Lorentz gamma factor.
double TwoShock_RarefactionSpeed | ( | double | w[], |
int | iside | ||
) |
Compute head or tail characteristic speeds enclosing the rarefaction fan:
[in] | w | a vector of primitive quantities containing the three velocities. |
[in] | iside(IN) | an integer specifying a left (-1) or right (+1) rarefaction wave. |
|
static |
Compute post shock quantities u1, dudp, zeta, for a given value of the post-shock pressure p1
Solve Riemann problem for the relativistic HD equations using the two-shock Riemann solver of Mignone et al. (2005).