PLUTO  4.4-patch2
Macros | Functions
hlld.c File Reference

HLLD Riemann solver for the relativistic MHD equations. More...

#include "pluto.h"

Macros

#define COUNT_FAILURES   NO
 

Functions

static double HLLD_Fstar (Riemann_State *, Riemann_State *, double)
 
static int HLLD_GetRiemannState (Riemann_State *, double, int)
 
static void HLLD_GetAState (Riemann_State *, double p)
 
static void HLLD_GetCState (Riemann_State *, Riemann_State *, double, double *)
 
static double HLLD_TotalPressure (double *)
 
void HLLD_Solver (const Sweep *sweep, int beg, int end, double *cmax, Grid *grid)
 

Detailed Description

Solve the Riemann problem for the relativistic MHD equations using the HLLD solver of Mignone, Ugliano & Bodo (2009). The solver requires the solution of a nonlinear equation (Eq. [48]) and the maximum number of iterations before convergence can be set using the macro MAX_ITER (default 20). A physical relevant solution is accepted if it satisfies the constraints by Eq. [54], implemented in the function HLLD_Fstar(). If this is not the case or if other failure conditions occur during the iteration loop, we switch to the HLL Riemann solver.

The macro COUNT_FAILURES can be set to YES if one wishes to count how many times the solver fails.

On input, it takes left and right primitive sweep vectors sweep->vL and sweep->vR 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:

Authors
A. Mignone (migno.nosp@m.ne@t.nosp@m.o.inf.nosp@m.n.it)
Date
Dec 03, 2019

Macro Definition Documentation

◆ COUNT_FAILURES

#define COUNT_FAILURES   NO

When set to YES, count number of failures and write the count to "hlld_fails.dat" (works in parallel as well). The number of failures is normalized to the total number of zones: nf/nz where nz is incremented during the main loop and nf is incremented when a failure occur.

Function Documentation

◆ HLLD_Fstar()

double HLLD_Fstar ( Riemann_State *  PaL,
Riemann_State *  PaR,
double  p 
)
static

Return the velocity difference across the contact mode as a function of the total pressure p.

◆ HLLD_GetAState()

void HLLD_GetAState ( Riemann_State *  Pa,
double  p 
)
static

Compute states aL and aR behind fast waves.

◆ HLLD_GetCState()

void HLLD_GetCState ( Riemann_State *  PaL,
Riemann_State *  PaR,
double  p,
double *  Uc 
)
static

Compute states cL and cR across contact mode.

◆ HLLD_GetRiemannState()

int HLLD_GetRiemannState ( Riemann_State *  Pv,
double  p,
int  side 
)
static

Express the sweep behind a wave as function of the total pressure p and the right hand side on the other side of the wave.

On output, return 1 if succesful, 0 if w < 0 is encountered. side = -1 : means left side = 1 : means right

When computing Bx, By and Bz, we use Eq. [21] with vx, vy and vz replaced by by Eq. [23,24,25]. The following short MAPLE script is used to verify the correctness.

restart;
A := R[mx] - lambda*R[E] + p*(1 - lambda^2):
G := R[By]*R[By] + R[Bz]*R[Bz]:
C := R[my]*R[By] + R[mz]*R[Bz]:
Q := - A - G + B[x]*B[x]*(1-lambda^2):
X := B[x]*(A*lambda*B[x]+ C) - (A + G)*(lambda*p + R[E]):
v[x] := (B[x]*(A*B[x] + lambda*C) - (A + G)*(p + R[mx]))/X:
v[y] := (Q*R[my] + R[By]*(C + B[x]*(lambda*R[mx] - R[E])))/X:
B[y] := (R[By] - B[x]*v[y])/(lambda - v[x]):
simplify(B[y]);

◆ HLLD_Solver()

void HLLD_Solver ( const Sweep sweep,
int  beg,
int  end,
double *  cmax,
Grid grid 
)

Solve the Riemann problem using the HLLD Riemann solver.

Parameters
[in,out]sweeppointer to Sweep structure
[in]beginitial grid index
[out]endfinal grid index
[out]cmax1D array of maximum characteristic speeds
[in]gridpointer to array of Grid structures.

◆ HLLD_TotalPressure()

double HLLD_TotalPressure ( double *  v)
static

Compute total pressure