Riemann problem for acoustics¶

This notebook illustrates some features of the exact solution to the 1-dimensional constant coefficient acoustics equations.

The notebook can be found in the clawpack apps respository:

$CLAW/apps/notebooks/riemann/Riemann_problem_acoustics.ipynb
See http://www.clawpack.org/apps.html.

In [1]:

%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np

The next cell imports a module containing a function that takes a Riemann problem (left state, right state, and approximate solver), and computes the Riemann solution, as well as functions to plot the solution in various forms.

In [3]:

#from clawpack.riemann import riemann_tools

# the version from clawpack.riemann needs updating, 
# for active development use the version in this repository instead:
import sys
sys.path.insert(0,'..')  
import riemann_tools

Acoustics: exact solution¶

We can use this to examine the exact solution of an acoustics Riemann problem.

This is a linear hyperbolic system of two equations for $q = [p, u]^T$, where $p$ is the pressure perturbation and $u$ is the velocity. The system is $q_t + Aq_x = 0$, where the coefficient matrix is

$$ A = \left[\begin{array}{cc}0&K\\1/\rho&0\end{array}\right], $$

where $\rho$ is the density and $K$ the bulk modulus. If we define the sound speed $c = \sqrt{K/\rho}$ and impedance $Z=\sqrt{K\rho}$, then the eigenvalues of the matrix are $s^1 = -c$ and $s^2 = +c$ and the corresponding eigenvectors are $$ r^1 = \left[\begin{array}{c}-Z\\1\end{array}\right], \qquad r^2 = \left[\begin{array}{c}Z\\1\end{array}\right]. $$

For arbitrary states $q_\ell$ and $q_r$, the Riemann solution consists of two waves propagating with velocities $\pm c$ with an intermediate state $q_m$ that is connected to $q_\ell$ by a multiple of $r^1$ and to $q_r$ by a multiple of $r^2$.

This Riemann solver can be solved by the PyClaw solver riemann.acoustics_1D_py.acoustics_1D:

In [4]:

from clawpack.riemann.acoustics_1D_py import acoustics_1D
solver = acoustics_1D

Set some problem data needed by the solver:¶

In [5]:

problem_data = {}

rho = 1.            # density
K = 4.              # bulk modulus
c = np.sqrt(K/rho)  # sound speed
Z = np.sqrt(K*rho)  # impedance

print("Density rho = %g, Bulk modulus K = %g" % (rho,K))
print("Sound speed = %g, Impedance = %g" % (c,Z))

problem_data['zz'] = Z 
problem_data['cc'] = c

Density rho = 1, Bulk modulus K = 4
Sound speed = 2, Impedance = 2

Set the left and right states, and solve the Riemann problem¶

In [6]:

q_l = np.array((1,4))  # Left state
q_r = np.array((3,7))  # Right state

states, s, riemann_eval = riemann_tools.riemann_solution(solver,q_l,q_r,\
                                            problem_data=problem_data, verbose=True)

States in Riemann solution:

$$\left [ \left[\begin{matrix}1.0\\4.0\end{matrix}\right], \quad \left[\begin{matrix}-1.0\\5.0\end{matrix}\right], \quad \left[\begin{matrix}3.0\\7.0\end{matrix}\right]\right ]$$

Waves (jumps between states):

$$\left [ \left[\begin{matrix}-2.0\\1.0\end{matrix}\right], \quad \left[\begin{matrix}4.0\\2.0\end{matrix}\right]\right ]$$

Speeds:

$$\left[\begin{matrix}-2.0 & 2.0\end{matrix}\right]$$

Fluctuations amdq, apdq:

$$\left [ \left[\begin{matrix}4.0 & -2.0\end{matrix}\right], \quad \left[\begin{matrix}8.0 & 4.0\end{matrix}\right]\right ]$$

Plot the states in the phase plane:¶

In [7]:

fig = plt.figure(figsize=(4,3))
ax = plt.axes()
riemann_tools.plot_phase(states, ax=ax, label_h='pressure', label_v='velocity')

Plot the waves in the x-t plane, and the solution at one particular time:¶

In [8]:

fig = riemann_tools.plot_riemann(states,s,riemann_eval,t=0.5)

Animate the Riemann solution:¶

In [9]:

riemann_tools.JSAnimate_plot_riemann(states,s,riemann_eval)

Out[9]: