Normal shock waves¶

%matplotlib inline
from matplotlib import pyplot as plt

import numpy as np

# Pint gives us some helpful unit conversion
from pint import UnitRegistry
ureg = UnitRegistry()
Q_ = ureg.Quantity # We will use this to construct quantities (value + unit)

# these lines are only for helping improve the display
import matplotlib_inline.backend_inline
matplotlib_inline.backend_inline.set_matplotlib_formats('pdf', 'png')
plt.rcParams['figure.dpi']= 300
plt.rcParams['savefig.dpi'] = 300
plt.rcParams['mathtext.fontset'] = 'cm'

Shock waves are finite pressure disturbances that correspond to a large change in properties over a short distance, on the order of a few free molecular paths of the gas. In contrast, sound waves are infinitesimal disturbances.

We can analyze the flow through a shock wave with a control volume, applying our existing conservation equations.

Fig. 1 Control volume around a shock.¶

Figure 1 shows a control volume around a shock wave, in a duct with varying area. The thickness of the control volume is very small, on the order of the thickness of the shock itself (so \(dx \sim 10^{-6}\) m). We will make the following assumptions about the flow:

steady, one-dimensional flow
adiabatic flow process: \(\delta q = 0\), and so \(ds_e = 0\)
no shaft work across the control volume: \(\delta w_s = 0\)
no potential change: \(dz = 0\)
constant area around the shock: \(A_1 = A_2\)

We can apply conservation of mass to this control volume:

\[ \rho_1 A_1 V_1 = \rho_2 A_2 V_2 \]

(23)¶\[ \rho_1 V_1 = \rho_2 V_2 \;, \]

conservation of energy:

\[\begin{gather*} h_{t1} + q = h_{t2} + w_s \\ h_{t1} = h_{t2} \end{gather*}\]

(24)¶\[ h_1 + \frac{V_1^2}{2} = h_2 + \frac{V_2^2}{2} \;, \]

and momentum:

\[\begin{gather*} \sum{F_x} = \dot{m} \left( V_{\text{out},x} - V_{\text{in}, x} \right) = \dot{m} \left( V_{2x} - V_{1x} \right) \\ p_1 A_1 - p_2 A_2 = ( p_1 - p_2 ) A = \dot{m} (V_2 - V_1) = \rho A V (V_2 - V_1) \\ p_1 - p_2 = (\rho V) V_2 - (\rho V) V_1 \\ p_1 - p_2 = \rho_2 V_2^2 - \rho_1 V_1^2 \end{gather*}\]

(25)¶\[ p_1 + \rho_1 V_1^2 = p_2 + \rho_2 V_2^2 \;. \]

For an arbitrary fluid, we have three equations: (23), (24), and (25). In a typical problem, we know the fluid and conditions before the shock, and want to find the conditions after the shock. Thus, our known variables are \(\rho_1\), \(p_1\), \(h_1\), and \(V_1\), while our unknown variables are \(\rho_2\), \(p_2\), \(h_2\), and \(V_2\). We need a fourth equation to close this system of equation: a property relation for the fluid, otherwise known as an equation of state.

Perfect gases¶

For ideal/perfect gases, we have the ideal gas equation of state, and can assume constant values for \(c_p\), \(c_v\), and \(\gamma\). We also have a convenient relationship for the speed of sound:

\[\begin{gather*} p = \rho R T \\ V = M a = M \sqrt{\gamma R T} \;. \end{gather*}\]

Incorporating these into Equation (23) we can obtain

(26)¶\[ \frac{p_1}{M_1}{\sqrt{T_1}} = \frac{p_2 M_2}{\sqrt{T_2}} \;. \]

With our stagnation relationship for temperature, \( T_t = T \left( 1 + (\gamma-1)M^2 / 2 \right)^2 \), Equation (24) becomes

(27)¶\[ T_1 \left( 1 + \frac{\gamma-1}{2} M_1^2 \right) = T_2 \left( 1 + \frac{\gamma-1}{2} M_2^2 \right) \;, \]

which applies both around a normal shock and also at any point in a general flow with no work or heat transfer. Lastly, incorporating our relationships for perfect gases into Equation (25) we obtain

(28)¶\[ p_1 \left( 1 + \gamma M_1^2 \right) = p_2 \left( 1 + \gamma M_2^2 \right) \;. \]

Now, we have three equations and three unknowns: \(M_2\), \(p_2\), and \(T_2\). 😎

However, it is a bit of a pain to solve this complicated system of equations every time. Fortunately, we can combine all three together and eliminate pressure and temperature completely!

\[ \left( \frac{1 + \gamma M_2^2}{1 + \gamma M_1^2} \right) \frac{M_1}{M_2} = \left( \frac{1 + \frac{\gamma-1}{2} M_2^2}{1 + \frac{\gamma-1}{2} M_1^2} \right)^{1/2} \;, \]

which can actually be solved to find \(M_2 = f(\gamma, M_1)\). A trivial solution to this equation is that \(M_1 = M_2\), where there is no shock. (What does this mean? Just that the governing equations we set up apply just fine for a flow where there is no shock—good news!)

For nontrivial solutions, with some painful algebra we can rearrange this equation into a recognizable form:

\[ A \left( M_2^2 \right)^2 + B M_2^2 + C = 0 \;, \]

where \(A, B, C = f(\gamma, M_1)\). This looks like a quadratic equation! The only physically viable solution to this equation is

(29)¶\[ M_2^2 = \frac{M_1^2 + \frac{2}{\gamma-1}}{\frac{2 \gamma}{\gamma-1} M_1^2 - 1} \;. \]

gamma = 1.4
mach1 = np.linspace(1.0, 4.0, num=50, endpoint=True)
mach2 = np.sqrt((mach1**2 + 2/(gamma-1))/(2*gamma*mach1**2/(gamma-1) - 1))

plt.plot(mach1, mach2)
plt.xlabel(r'$M_1$')
plt.ylabel(r'$M_2$')
plt.title('Downstream Mach number versus upstream Mach number for normal shocks')
plt.grid(True)
plt.tight_layout()
plt.show()

A typical problem flow is that we have \(\gamma\) and \(M_1\), use these to determine \(M_2\), and then use Equations (27) and (28) to obtain \(T_2\) and \(p_2\). More convenient versions of these equations are

(30)¶\[ \frac{T_2}{T_1} = \frac{1 + \frac{\gamma-1}{2} M_1^2}{1 + \frac{\gamma-1}{2} M_2^2} \]

(31)¶\[ \frac{p_2}{p_1} = \frac{1 + \gamma M_1^2}{1 + \gamma M_2^2} = \frac{2 \gamma M_1^2 - \gamma + 1}{\gamma + 1} \]

(32)¶\[ \frac{\rho_2}{\rho_1} = \frac{(\gamma+1) M_1^2}{(\gamma-1) M_1^2 + 2} \]

Note

Since \(M_1 > 1\), \(M_2 < 1\) always for a normal shock.

Gas Dynamics notes

Normal shock waves¶

Perfect gases¶