Shallow water equations with a tracer

We can augment the shallow water equations discussed in Shallow_water.html with a quantity $\phi(x,t)$ that measures the concentration of a tracer that is advected with the fluid motion (but has no influence on the fluid dynamics). If $\phi$ is measured in units of mass per unit volume (which is really per unit area of the vertical cross section in this one-dimensional example), then the mass per unit length along the $x$ axis is given by $h(x,t)\phi(x,t)$. The quantity $\phi$ satisfies the variable coefficient advection equation in advective (nonconservative) form: $$ \phi_t(x,t) + u(x,t)\phi_x(x,t) = 0. $$ This equation is also called the "color equation", since we can think of $\phi$ as measuring the concentration of a dye that changes the color of the water. We will use this interpretation in the plots below. By setting the intitial conditions $\phi(x,0)$ to be piecewise constant with different values corresponding to different colors, we can visualize the motion of the fluid better. We will use a dark shade of blue for the water that is initially to the left of a dam at $x=0$ and a lighter shade of blue for the water that is initially to the right.

The quantity $h\phi$ satisfies the conservative form of the advection equation, $$ (h\phi)_t + (uh\phi)_x = 0. $$ This can be derived by combining the color equation with the conservation of mass equation $h_t +(hu)_x = 0$. Since $h\phi$ measures the concentration of dye per unit length in $x$, and we assume molecules of dye are not created or destroyed, it makes sense that this is the conserved quantity.

The full system of equations in conservation form, \begin{align} \label{SW-tracer} q_t(x,t) + f(q(x,t))_x = 0, \end{align} thus takes the form \begin{align*} \begin{split} h_t + (hu)_x &=0,\\ (hu)_t + \left(hu^2 + \frac 1 2 gh^2\right)_x &=0,\\ (h\phi)_t + (uh\phi)_x &= 0. \end{split} \end{align*}

The Riemann solution for this system has 3 waves. Two of these waves, with speeds $u\pm c$ (where $c=\sqrt{gh}$) are just the same waves we studied in the shallow water system. We refer to these as the 1-wave and 2-wave, to maintain consistency with the last chapter. The additional wave appearing in this system (the tracer wave) moves at the fluid velocity $u$ and carries only a jump in the tracer concentration $\phi$. It is our first example of what is known as a contact discontinuity or linearly degenerate wave. Because the fluid velocity is constant across this wave, tracer characteristics travel parallel to the wave on either side. Thus this family does not admit shock or rarefaction waves. We will see another example of linear degeneracy in Euler and give a more general definition at that point.

The 1- and 2-waves can of course each be either a shock wave or rarefaction wave depending on the initial data, as explained in Shallow_water.html, since the tracer has no effect on the structure of these waves. Below we consider a dam break problem in which case there is one of each.

In [1]:
%matplotlib inline
In [2]:
%config InlineBackend.figure_format = 'svg'
from ipywidgets import widgets, fixed
from utils.snapshot_widgets import interact

from exact_solvers import shallow_water
demo_plot = shallow_water.make_demo_plot_function
Will create static figures with single value of parameters

If you wish to examine the Python code for this chapter, see

Dam break problem

Here we plot the solution to the dam break Riemann problem for the shallow water equations with the addition of a tracer. The tracer is now also taken to be piecewise constant, with light blue representing the value $\phi_l$ taken in the water the left of the dam, while dark blue represents the value $\phi_r$ in the water on the right side of the dam.

The structure of the depth and velocity are exactly the same as seen and discussed in detail in Shallow_water.html, and the value of $\phi$ is constant across the 1-rarefaction and 2-shock waves. The discontinuity in $\phi$ shows the propagation of the contact wave in the Riemann solution, the contact discontinuity across which $\phi$ is discontinuous while $h$ and $u$ are continuous.

In [3]:
shallow_water.plot_riemann_SW(h_l=3, h_r=1, u_l=0., u_r=0., tracer=True)

In the live notebook you can select which set of characteristics to plot. Note that the characteristics corresponding to the tracer are equivalent to the particle paths of the water that were shown in similar examples in Shallow_water.html. This is because the characteristic speed is equal to $u$, the fluid velocity.

Hyperbolic structure

The canonical form (\ref{SW-tracer}) arises by defining \begin{align} q & = \begin{pmatrix} h \\ hu \\ h \phi \end{pmatrix}, & f & = \begin{pmatrix} hu \\ hu^2 + \frac{1}{2}gh^2 \\ uh\phi\end{pmatrix}. \end{align}
In terms of the conserved quantities, the flux is
\begin{align} f(q) & = \begin{pmatrix} q_2 \\ q_2^2/q_1 + \frac{1}{2}g q_1^2 \\ q_3 q_2/q_1 \end{pmatrix}. \end{align}
Thus the flux Jacobian is
\begin{align} f'(q) & = \begin{pmatrix} 0 & 1 & 0 \\ -(q_2/q_1)^2 + g q_1 & 2 q_2/q_1 & 0 \\ -q_3 q_2/q_1^2 & q_3/q_1 & q_2/q_1 \end{pmatrix} = \begin{pmatrix} 0 & 1 & 0 \\ -u^2 + g h & 2 u & 0 \\ -u\phi & \phi & u \end{pmatrix}. \end{align}
Its eigenvalues are
\begin{align*} \lambda_1 & = u - \sqrt{gh} & \lambda_\text{tracer} & = u & \lambda_2 & = u + \sqrt{gh}, \end{align*}
with corresponding eigenvectors
\begin{align*} r_1 & = \begin{bmatrix} 1 \\ u-\sqrt{gh} \\ \phi \end{bmatrix} & r_\text{tracer} & = \begin{bmatrix} 0 \\ 0 \\ 1 \end{bmatrix} & r_2 & = \begin{bmatrix} 1 \\ u+\sqrt{gh} \\ \phi \end{bmatrix} \end{align*}
We see immediately that the tracer wave carries only a jump in the tracer itself. It may seem puzzling that the third entry of $r_1$ and $r_3$ is nonzero. Notice that each of these waves carries a jump in $h$; since $\phi$ is constant across each of these waves, there must be a corresponding jump in $q_3=\phi h$.

Linear degeneracy

The speed of the tracer wave depends only on $u$, and $u$ is constant across the tracer wave. This means that tracer characteristics can neither converge nor diverge due to a jump in $\phi$, but rather run parallel to it. Characteristic fields with this property are referred to as being linearly degenerate since the Riemann solution to the linear advection equation has this same property, that the characteristic speed is independent of changes in the variable that jumps across the wave.

If we compute the gradient of the eigenvalue $\lambda_\text{tracer}$ by differentiating $u = hu/h$ with respect to each of the component of $q=(h, hu, h\phi)$, we obtain \begin{align} \label{SW:gradu} \nabla \lambda_\text{tracer} = \begin{bmatrix} -hu/h^2\\ 1/h \\ 0\end{bmatrix}. \end{align} This vector is orthogonal to the corresponding eigenvector, \begin{align} \label{SW:lindeg} \nabla \lambda_\text{tracer} \cdot r_\text{tracer} \equiv 0 \end{align} for all values of $q$. This is the general mathematical condition for a linearly degenerate field, since it shows that the eigenvalue is unchanged as we move along integral curves of the corresponding eigenvalue.

We will see another example of a linearly degenerate characteristic family in the next chapter, Euler.html, where we will see that the Riemann solution to the Euler equations has a very similar structure to that of the shallow water equations with a tracer. In general it consists of two nonlinear acoustic waves that can be either shock or rarefaction waves, depending on the initial data, and an intermediate contact discontinuity that propagates with the intermediate fluid velocity, and across which the fluid velocity is constant.