Riemann Solver Package

This package contains all of the Python-based Riemann solvers. Each module solves the Riemann solver for a particular system of hyperbolic equations. The solvers all have a common function signature:

rp_<name>_<dim>d(q_l,q_r,aux_l,aux_r,problem_data)

with <name> replaced with the appropriate solver name and <dim> with the appropriate dimension.

Input:
  • q_l - (ndarray(...,num_eqn)) Contains the left states of the Riemann problem
  • q_r - (ndarray(...,num_eqn)) Contains the right states of the Riemann problem
  • aux_l - (ndarray(...,num_aux)) Contains the left values of the auxiliary array
  • aux_r - (ndarray(...,num_aux)) Contains the right values oft he auxiliary array
  • problem_data - (dict) Dictionary containing miscellaneous data which is
    usually problem dependent.
Output:
  • wave - (ndarray(...,num_eqn,num_waves)) Contains the resulting waves from the cell
    edge
  • s - (ndarray(...,num_waves)) Speeds of each wave
  • amdq - (ndarray(...,num_eqn)) Left going fluctuation
  • apdq - (ndarray(...,num_eqn)) Right going fluctuation

Except for problem_data, all of the input and output values are arrays whose elements represent grid values with locations indicated by the following scheme

Indexing works like this:  here num_ghost=2 as an example
 0     1     2     3     4     mx+num_ghost-2     mx+num_ghost      mx+num_ghost+2
             |                        mx+num_ghost-1 |  mx+num_ghost+1
 |     |     |     |     |   ...   |     |     |     |     |
    0     1  |  2     3            mx+num_ghost-2    |mx+num_ghost
                                          mx+num_ghost-1   mx+num_ghost+1

The top indices represent the values that are located on the grid
cell boundaries such as waves, s and other Riemann problem values,
the bottom for the cell centered values.  In particular the ith grid cell
boundary has the following related information:
                  i-1         i         i+1
                   |          |          |
                   |   i-1    |     i    |
                   |          |          |
Again, grid cell boundary quantities are at the top, cell centered
values are in the cell.

Note

The values q_l[i], q_r[i] are the left and right states, respectively, of the ith Riemann problem. This convention is different than that used in the Fortran Riemann solvers, where q_l[i], q_r[i] are the values at the left and right edges of a cell.

All of the return values (waves, speeds, and fluctuations) are indexed by cell edge (Riemann problem being solved), with s[i] referring to the wave speed at interface $i-1/2$. This follows the same convention used in the Fortran solvers.

See [LeVeque_book_2002] for more details.

List of available Riemann solvers:

Acoustics

Riemann solvers for constant coefficient acoustics.

\[q_t + A q_x = 0\]

where

\[\begin{split}q(x,t) = \left [ \begin{array}{c} p(x,t) \\ u(x,t) \end{array} \right ]\end{split}\]

and the coefficient matrix is

\[\begin{split}A = \left [\begin{matrix} 0 & K\\ 1/\rho & 0 \end{matrix} \right ]\end{split}\]

The parameters \(\rho =\) density and \(K =\) bulk modulus are used to calculate the impedence \(= Z\) and speed of sound = c.

Authors:Kyle T. Mandli (2009-02-03): Initial version
clawpack.riemann.acoustics_1D_py.acoustics_1D(q_l, q_r, aux_l, aux_r, problem_data)

Basic 1d acoustics riemann solver, with interleaved arrays

problem_data is expected to contain -
  • zz - (float) Impedence
  • cc - (float) Speed of sound

See Riemann Solver Package for more details.

Version:1.0 (2009-02-03)

Advection

Simple advection Riemann solvers

Basic advection Riemann solvers of the form (1d)

\[q_t + A q_x = 0.\]
Authors:Kyle T. Mandli (2008-2-20): Initial version
clawpack.riemann.advection_1D_py.advection_1D(q_l, q_r, aux_l, aux_r, problem_data)

Basic 1d advection riemann solver

problem_data should contain -
  • u - (float) Determines advection speed

See Riemann Solver Package for more details.

Version:1.0 (2008-2-20)

Burgers Equation

Riemann solvers for Burgers equation

\[u_t + \left ( \frac{1}{2} u^2 \right)_x = 0\]
Authors:Kyle T. Mandli (2009-2-4): Initial version
clawpack.riemann.burgers_1D_py.burgers_1D(q_l, q_r, aux_l, aux_r, problem_data)

Riemann solver for Burgers equation in 1d

problem_data should contain -
  • efix - (bool) Whether a entropy fix should be used, if not present, false is assumed

See Riemann Solver Package for more details.

Version:1.0 (2009-2-4)

Euler Equations

Riemann solvers for the Euler equations

This module contains Riemann solvers for the Euler equations which have the form (in 1d):

\[q_t + f(q)_x = 0\]

where

\[\begin{split}q(x,t) = \left [ \begin{array}{c} \rho \\ \rho u \\ E \end{array} \right ],\end{split}\]

the flux function is

\[\begin{split}f(q) = \left [ \begin{array}{c} \rho u \\ \rho u^2 + p \\ u(E+p) \end{array}\right ].\end{split}\]

and \(\rho\) is the density, \(u\) the velocity, \(E\) is the energy and \(p\) is the pressure.

Unless otherwise noted, the ideal gas equation of state is used:

\[E = (\gamma - 1) \left (E - \frac{1}{2}\rho u^2 \right)\]
Authors:Kyle T. Mandli (2009-06-26): Initial version Kyle T. Mandli (2011-03-28): Interleaved version
clawpack.riemann.euler_1D_py.euler_exact_1D(q_l, q_r, aux_l, aux_r, problem_data)

Exact euler Riemann solver

Warning

This solver has not been implemented.

clawpack.riemann.euler_1D_py.euler_hll_1D(q_l, q_r, aux_l, aux_r, problem_data)

HLL euler solver

W_1 = Q_hat - Q_l    s_1 = min(u_l-c_l,u_l+c_l,lambda_roe_1,lambda_roe_2)
W_2 = Q_r - Q_hat    s_2 = max(u_r-c_r,u_r+c_r,lambda_roe_1,lambda_roe_2)

Q_hat = ( f(q_r) - f(q_l) - s_2 * q_r + s_1 * q_l ) / (s_1 - s_2)
problem_data should contain:
  • gamma - (float) Ratio of the heat capacities
  • gamma1 - (float) \(1 - \gamma\)
Version:1.0 (2014-03-04)
clawpack.riemann.euler_1D_py.euler_hllc_1D(q_l, q_r, aux_l, aux_r, problem_data)

HLLC Euler solver

System Message: WARNING/2 (/Users/rjl/git/clawpack/clawpack/riemann/euler_1D_py.py:docstring of clawpack.riemann.euler_1D_py.euler_hllc_1D, line 3)

Literal block expected; none found.

W_1 = q_hat_l - q_l s_1 = min(u_l-c_l,u_l+c_l,lambda_roe_1,lambda_roe_2) W_2 = q_hat_r - q_hat_l s_2 = s_m W_3 = q_r - q_hat_r s_3 = max(u_r-c_r,u_r+c_r,lambda_roe_1,lambda_roe_2)

s_m = (p_r - p_l + rho_l*u_l*(s_l - u_l) - rho_r*u_r*(s_r - u_r))
/ (rho_l*(s_l-u_l) - rho_r*(s_r - u_r))

# left middle state q_hat_l[0,:] = rho_l*(s_l - u_l)/(s_l - s_m) q_hat_l[1,:] = rho_l*(s_l - u_l)/(s_l - s_m)*s_m q_hat_l[2,:] = rho_l*(s_l - u_l)/(s_l - s_m)

System Message: ERROR/3 (/Users/rjl/git/clawpack/clawpack/riemann/euler_1D_py.py:docstring of clawpack.riemann.euler_1D_py.euler_hllc_1D, line 14)

Unexpected indentation.

*(E_l/rho_l + (s_m - u_l)*(s_m + p_l/(rho_l*(s_l - u_l))))

System Message: WARNING/2 (/Users/rjl/git/clawpack/clawpack/riemann/euler_1D_py.py:docstring of clawpack.riemann.euler_1D_py.euler_hllc_1D, line 14); backlink

Inline emphasis start-string without end-string.

# right middle state q_hat_r[0,:] = rho_r*(s_r - u_r)/(s_r - s_m) q_hat_r[1,:] = rho_r*(s_r - u_r)/(s_r - s_m)*s_m q_hat_r[2,:] = rho_r*(s_r - u_r)/(s_r - s_m)

System Message: ERROR/3 (/Users/rjl/git/clawpack/clawpack/riemann/euler_1D_py.py:docstring of clawpack.riemann.euler_1D_py.euler_hllc_1D, line 20)

Unexpected indentation.

*(E_r/rho_r + (s_m - u_r)*(s_m + p_r/(rho_r*(s_r - u_r))))

System Message: WARNING/2 (/Users/rjl/git/clawpack/clawpack/riemann/euler_1D_py.py:docstring of clawpack.riemann.euler_1D_py.euler_hllc_1D, line 20); backlink

Inline emphasis start-string without end-string.

problem_data should contain: - gamma - (float) Ratio of specific heat capacities - gamma1 - (float) \(\gamma - 1\)

:Version 1.0 (2015-11-18)

clawpack.riemann.euler_1D_py.euler_roe_1D(q_l, q_r, aux_l, aux_r, problem_data)

Roe Euler solver in 1d

aug_global should contain -
  • gamma - (float) Ratio of the heat capacities
  • gamma1 - (float) \(1 - \gamma\)
  • efix - (bool) Whether to use an entropy fix or not

See Riemann Solver Package for more details.

Version:1.0 (2009-6-26)

Shallow Water Equations

Riemann solvers for the shallow water equations.

The available solvers are:
  • Roe - Use Roe averages to caluclate the solution to the Riemann problem
  • HLL - Use a HLL solver
  • Exact - Use a newton iteration to calculate the exact solution to the
    Riemann problem
\[q_t + f(q)_x = 0\]

where

\[\begin{split}q(x,t) = \left [ \begin{array}{c} h \\ h u \end{array} \right ],\end{split}\]

the flux function is

\[\begin{split}f(q) = \left [ \begin{array}{c} h u \\ hu^2 + 1/2 g h^2 \end{array}\right ].\end{split}\]

and \(h\) is the water column height, \(u\) the velocity and \(g\) is the gravitational acceleration.

clawpack.riemann.shallow_1D_py.shallow_exact_1D(q_l, q_r, aux_l, aux_r, problem_data)

Exact shallow water Riemann solver

Warning

This solver has not been implemented.

clawpack.riemann.shallow_1D_py.shallow_fwave_1d(q_l, q_r, aux_l, aux_r, problem_data)

Shallow water Riemann solver using fwaves

Also includes support for bathymetry but be wary if you think you might have dry states as this has not been tested.

problem_data should contain:
  • grav - (float) Gravitational constant
  • dry_tolerance - (float) Set velocities to zero if h is below this tolerance.
  • sea_level - (float) Datum from which the dry-state is calculated.
Version:1.0 (2014-09-05)
Version:2.0 (2017-03-07)
clawpack.riemann.shallow_1D_py.shallow_hll_1D(q_l, q_r, aux_l, aux_r, problem_data)

HLL shallow water solver

W_1 = Q_hat - Q_l    s_1 = min(u_l-c_l,u_l+c_l,lambda_roe_1,lambda_roe_2)
W_2 = Q_r - Q_hat    s_2 = max(u_r-c_r,u_r+c_r,lambda_roe_1,lambda_roe_2)

Q_hat = ( f(q_r) - f(q_l) - s_2 * q_r + s_1 * q_l ) / (s_1 - s_2)
problem_data should contain:
  • g - (float) Gravitational constant
Version:1.0 (2009-02-05)
clawpack.riemann.shallow_1D_py.shallow_roe_1D(q_l, q_r, aux_l, aux_r, problem_data)

Roe shallow water solver in 1d:

ubar = (sqrt(u_l) + sqrt(u_r)) / (sqrt(h_l) + sqrt(h_r))
cbar = sqrt( 0.5 * g * (h_l + h_r))

W_1 = |      1      |  s_1 = ubar - cbar
      | ubar - cbar |

W_2 = |      1      |  s_1 = ubar + cbar
      | ubar + cbar |

a1 = 0.5 * ( - delta_hu + (ubar + cbar) * delta_h ) / cbar
a2 = 0.5 * (   delta_hu - (ubar - cbar) * delta_h ) / cbar
problem_data should contain:
  • g - (float) Gravitational constant
  • efix - (bool) Boolean as to whether a entropy fix should be used, if not present, false is assumed
Version:1.0 (2009-02-05)