Modeling Shock Waves

Sultan Sial

The modeling of shock waves illustrates themes that are common in the field of Applied Mathematics:

Pictures for the impatient

We will illustrate these themes by looking at a simple case - the Sod shock tube. Click here to see the mass density in the tube after the shocks have developed. Click here to see the velocities in the shock tube after the shocks have developed.

Probabilistic Viewpoint

There are a number of approaches to thinking about and modeling a fluid like air or water.

One way is to try and describe the motion of each individual particle of the fluid. This leads to the field of molecular dynamics. A problem with this method is that it is very computationally intensive to track even a few thousand molecules let alone the millions in even a very small volume.

Another approach is to assign a probability density function f to the fluid that describes fluid.

In statistical mechanics such a function takes a form f(u,v,w) where u, v, w are the x-, y-, and z-velocities of individual molecules. f tells us what fraction of the molecules of the fluid have the velocities u, v and w.

What form should the function f take? This questions leads to our first assumption.

Maxwell-Boltzmann Distribution

It is a principle of statistical mechanics that a system will seek to maximize 'disorder' or entropy. For an ideal gas, the distribution that maximizes the entropy is the Maxwell-Boltzmann distribution. So we could assume that, at any point in fluid, the fluid can be described by the Maxwell-Boltzmann distribution. This is of course false. However, we have to start somewhere and the Maxwell-Boltzmann distribution may be a good approximation in some cases.

Finite-Volume Method

There are several methods for describing the evolution of a fluid. The Finite-Volume methods take the volume of interest and divide it up into a number of smaller volumes.

We make the assumption that within each small finite volume the fluid will have the same mass, momentum, and energy. Again this is false, but with enough finite volumes is often a good approximation.

Numerical Interpolations

We can specify the mass, momentum and energy in each cell at the beginning of the simulation. But what about the values at the interfaces of the cells? At this point some sort of numerical interpolation is required. Numerical interpolation techniques attempt to make estimates for the value of a function when the exact value is difficult, expensive, or impossible to determine. Using the values in neighbouring finite volumes one can estimate the values at th interfaces between the finite volumes.


The rate at which any quantity s is flowing in the x-direction is given by integrating the function usf over all values of from minus infinity to infinity. The fluxes in the other two directions are given by integrating vsf and wsf.

We require one more false assumption. The function f, which we assume is the Maxwell-Boltzmann distribution, is a function of the average momentum at the point of interest. Assume that over a small time step T the average momentum does not change. Then, the amount of mass, momentum or energy transferred between two finite volumes will be given by the flux of that quantity times T.


With these assumptions it is possible to watch a fluid evolve. We begin a simulation by specifying the mass, momentum and energy in each finite volume. We then interpolate the values of these quantities at the interfaces between the finite volumes. These values give us the Maxwell-Boltzmann distribution at the interfaces from which the fluxes can be calculated. The fluxes tell us how much mass, momentum, and energy is transferred between the finite volumes in each time step. At the end of the time step we can add these fluxes to each finite volume to get new values for the mass, momentum, and energy. We then repeat the entire process over and over again to see how the fluid evolves.


An illustration of a finite volume technique is the simulation of the Sod shock tube. This is a standard test case for which analytical techniques give an answer against which numerical techniques can be tested.

The Sod test case consists of a one dimensional tube in which one half of the tube has different initial conditions than the other half. initially, the left half of the tube has a mass density of 1 and an energy of 2.5. The right side of tube initially has a mass density of 0.125 and an energy of 0.25.

As the flow develops shock waves form.

The results shown below are from a FORTRAN90 program run on an IBM RISC6000 workstation. Since the code is essentially parallel it runs well on super-computers also.

Initially, the shock tube has a high and a low density region. As the flow develops, regions of intermediate density develop as well. This is shown in the digram below.

Initially, the fluid is at rest. Then a flow develops in the middle of the tube which moves towards the right.