Conservation Laws and Shocks

6B.1 First-Order Equations and Characteristics ended at the moment Burgers' equation breaks: characteristics carrying different values cross, and no smooth solution continues. Physics does not stop there. A gas compressed by a piston forms a shock wave; cars arriving at the back of a queue form a sharp front where speed drops abruptly; a flood wave steepens into a bore. In each case the quantity that is conserved (mass, cars, water) is still conserved after the jump forms. The mathematical question is how to continue the solution past the breaking time in a way that keeps conservation, and the answer has three parts.

First, write the equation in integral form, which makes sense for discontinuous functions. That determines the speed of a jump: the Rankine–Hugoniot condition. Second, the integral form alone allows too many solutions, some of them unphysical, so an extra entropy condition is imposed to select the right one. Third, the selected solution is the limit of solutions of the equation with a tiny amount of viscosity added, which is the physical justification for the entropy condition, and for Burgers' equation this limit can be computed exactly, because the viscous equation turns into the heat equation.

The theme of the chapter, that a weak formulation plus a selection principle continues a solution through a singularity in a unique way, is the same theme as viscosity solutions of Hamilton–Jacobi equations (6C.9 Hamilton–Jacobi Equations and Viscosity Solutions), and, at a much higher level, as the weak Ricci flows that continue through singularities (12C.9 After Perelman).

By the end of this chapter you will be able to:

• define integral (weak) solutions of a scalar conservation law, and derive the Rankine–Hugoniot condition for the speed of a shock;
• show by example that integral solutions are not unique, and state and use the Lax entropy condition;
• solve the Riemann problem for a convex flux in full: shocks and rarefaction waves;
• turn viscous Burgers into the heat equation by the Hopf–Cole transformation, and compute the viscous shock profile and its width;
• compute shock speeds in a traffic-flow model, and the entropy dissipated by a shock.

2.1Traffic jams as shock waves

2.2Integral solutions

Consider a scalar conservation law in one space dimension,

with $F$ smooth (the flux). It is the differential form of the statement that for every interval $[a, b]$,

the amount of $u$ in $[a, b]$ changes only by what flows in at $a$ and out at $b$ (6A.1 What a PDE Is). Burgers' equation is the case $F(u) = \frac12u^2$, since $(\frac12u^2)_x = uu_x$.

To allow discontinuous $u$, multiply the equation by a smooth test function $v$ with compact support in $\mathbb{R}\times[0, \infty)$ (it need not vanish at $t = 0$), integrate, and move all derivatives onto $v$. For a smooth solution, integrating by parts in $t$ and in $x$,

The right side makes sense for any bounded measurable $u$.

A $C^1$ function is an integral solution exactly when it is a classical solution with $u(\cdot, 0) = g$ (Exercise 2.10). The interesting integral solutions are the ones with jumps.

2.2.1The Rankine–Hugoniot condition

Suppose $u$ is smooth on either side of a curve $x = s(t)$ in space-time and jumps across it, with limits $u_l$ from the left and $u_r$ from the right. What does the integral formulation say about the curve?

For Burgers' equation $F(u) = \frac12u^2$, the jump condition gives $\sigma(u_l - u_r) = \frac12(u_l^2 - u_r^2)$, so

a Burgers shock moves at the average of the values on its two sides. For the traffic flux it gives the speed of the back of a queue (Exercise 2.14).

2.3Non-uniqueness and the entropy condition

The integral formulation allows jumps, and it allows too many.

Which is right? Look at the characteristics. Characteristics of Burgers' equation travel at speed $u$, so on the left (where $u = 0$) they are vertical lines, and on the right (where $u = 1$) lines of slope $1$. In $u_1$, characteristics emerge from the shock on both sides: the value $u$ along them does not come from the initial data but is created at the shock. That violates causality: information should flow forward from the data, and into a shock, not out of it. In $u_2$, every point is reached by a characteristic from the initial line, the fan filling the wedge $0 < x < t$ with characteristics issuing from the origin.

For a convex flux, $F'$ is increasing, and by the mean value theorem $\sigma = \frac{F(u_l) - F(u_r)}{u_l - u_r} = F'(\xi)$ for some $\xi$ between $u_l$ and $u_r$. So $F'(u_l) > \sigma > F'(u_r)$ holds exactly when $u_l > u_r$. In Example 2.3 the shock has $u_l = 0 < 1 = u_r$, so it violates the entropy condition and is rejected. The rarefaction wave is the physically correct solution.

For general bounded solutions, which may have more complicated singularities than single shocks, Olga Oleinik's version of the condition (1957) is a one-sided bound on how fast $u$ can increase:

for a constant $C$. An upward jump violates it; a downward jump does not. A bounded integral solution satisfying this condition is called an entropy solution.

Existence is Evans, section 3.4.2, by way of the Hopf–Lax formula for the Hamilton–Jacobi equation $w_t + F(w_x) = 0$, whose derivative $u = w_x$ solves the conservation law (6C.9 Hamilton–Jacobi Equations and Viscosity Solutions proves Hopf–Lax). Uniqueness is Evans, section 3.4.3 (Oleinik's theorem), proved by a duality argument with a backward transport equation. For general fluxes, and in several space dimensions, Stanislav Kruzhkov (1970) gave the definitive theory, with the entropy condition in the form of inequalities $\eta(u)_t + q(u)_x \leq 0$ for every convex entropy $\eta$ with entropy flux $q' = \eta'F'$ (Exercise 2.16).

2.4The Riemann problem

The Riemann problem has the simplest discontinuous data, a single jump:

Its solution is the building block of the theory and of numerical methods: Godunov's scheme (1959) solves a Riemann problem at every cell boundary at every time step.

The rarefaction wave is self-similar: it depends only on $x/t$, as the data are unchanged by the scaling $(x, t) \mapsto (\lambda x, \lambda t)$ (6A.1 What a PDE Is made the same argument for parabolic scaling). It is the fan of characteristics issuing from the jump, each carrying a value between $u_l$ and $u_r$.

2.5Viscosity: the Hopf–Cole transformation

Real fluids have viscosity, and real drivers look ahead. Adding a small diffusion term to the conservation law,

gives a parabolic equation, which has smooth solutions for all time (6A.8 Nonlinear Parabolic Equations). The entropy solution should be the limit as $\varepsilon \to 0$: the vanishing viscosity method. For Burgers' equation this can be done exactly.

Since the heat equation can be solved explicitly (6A.3 The Heat Equation on ℝⁿ), so can viscous Burgers: with $\varphi(x, 0) = \exp\big(-\frac{1}{2\varepsilon}\int_0^xg\big)$,

(Exercise 2.13). As $\varepsilon \to 0$ the weights concentrate at the minimiser of the bracket, by Laplace's method (the largest value of $e^{-K/\varepsilon}$ dominates), and $u \to \frac{x - y(x, t)}{t}$: exactly the Lax–Oleinik formula for Burgers, where $L(q) = \frac12q^2$. The vanishing-viscosity limit is the entropy solution. (Evans, section 4.5.2, gives Laplace's method; the conclusion holds where the minimiser is unique.)

2.5.1The viscous shock profile

What does a shock look like with a little viscosity? Look for a travelling wave $u(x, t) = U(x - \sigma t)$ connecting $u_l$ at $-\infty$ to $u_r$ at $+\infty$.

The shock has width of order $\frac{\varepsilon}{u_l - u_r}$: as viscosity vanishes the profile sharpens into a jump, and stronger shocks are thinner. In air, shock waves are only a few mean free paths of the gas molecules thick, a fraction of a micrometre at ordinary pressures, which is why they look like true discontinuities. And the entropy condition has a physical meaning: an admissible shock is one that is the limit of a viscous profile.

2.6Long-time behaviour

Entropy solutions forget their data in a precise way. For Burgers' equation with $g$ integrable and compactly supported, the solution approaches, as $t \to \infty$, an N-wave: a function equal to $\frac xt$ on an interval $(-\sqrt{2pt}, \sqrt{2qt})$ and $0$ outside, with jumps at both ends, where $p$ and $q$ are determined by the data (the largest negative and positive partial integrals of $g$). Its size decays like $t^{-1/2}$ (Evans, section 3.4.5). The decay comes from shocks: each shock dissipates entropy at a rate cubic in its strength (Exercise 2.16), so strong shocks weaken quickly. This is the signature of sonic booms: the pressure signature of a supersonic aircraft, far from it, has the shape of an N.

2.7Systems of conservation laws

Gas dynamics, shallow water and traffic models with two variables are systems $\mathbf u_t + \mathbf F(\mathbf u)_x = 0$, with $\mathbf u(x, t) \in \mathbb{R}^m$. The system is hyperbolic if the matrix $D\mathbf F(\mathbf u)$ has real eigenvalues $\lambda_1(\mathbf u) < \dots < \lambda_m(\mathbf u)$ with a basis of eigenvectors: then a small disturbance splits into $m$ waves travelling at the speeds $\lambda_k$, the characteristic speeds. For one-dimensional gas dynamics they are $v - c$, $v$ and $v + c$, with $v$ the fluid velocity and $c$ the speed of sound: sound waves travel at $\pm c$ relative to the fluid, and the entropy is carried with it.

The Rankine–Hugoniot condition becomes the vector equation $\sigma[\![\mathbf u]\!] = [\![\mathbf F(\mathbf u)]\!]$, which says $[\![\mathbf u]\!]$ is an eigenvector of an averaged Jacobian; the Lax condition requires that exactly the right number of characteristics run into the shock. Peter Lax showed in 1957 that the Riemann problem for a strictly hyperbolic system, with each characteristic field either genuinely nonlinear or linearly degenerate, is uniquely solvable for small jumps by a composition of shocks, rarefactions and contact discontinuities (Evans, section 11.2). For two equations, Riemann invariants, functions constant along one family of characteristics, reduce the system to a pair of transport equations, and this was how Riemann analysed gas dynamics in 1860 (Evans, section 11.3). Global existence for systems with small total variation is a deep theorem of James Glimm (1965), and uniqueness of such solutions was proved by Alberto Bressan and collaborators around 2000.

2.8History

Bernhard Riemann's 1860 paper on the propagation of sound waves of finite amplitude introduced the Riemann problem and Riemann invariants. William Rankine (1870) and Pierre-Henri Hugoniot (1887–89) derived the jump conditions. Eberhard Hopf (1950) and Julian Cole (1951) independently found the transformation of viscous Burgers into the heat equation. Peter Lax's entropy condition and his solution of the Riemann problem for systems date from 1957, and Olga Oleinik's uniqueness theorem for scalar laws from 1957 as well. Sergei Godunov's scheme appeared in 1959, James Glimm's existence theorem in 1965, and Stanislav Kruzhkov's theory of entropy solutions in several dimensions in 1970. Lighthill and Whitham's traffic model was published in 1955 and Richards' in 1956.

2.9Exercises

Solution

Integrating by parts in the definition (allowed for $C^1$ $u$), the integral condition becomes $\iint(u_t + F(u)_x)v + \int(u(x, 0) - g(x))v(x, 0)\,dx = 0$ for all $v$. Test functions supported in $t > 0$ give $u_t + F(u)_x = 0$ (a continuous function orthogonal to all test functions vanishes); then arbitrary $v(\cdot, 0)$ gives $u(\cdot, 0) = g$. The converse is the same computation read backwards.

Solution

$u_l = 1 > 0 = u_r$, so the solution is a shock at speed $\frac12$: $u = 1$ for $x < \frac t2$, $u = 0$ for $x > \frac t2$. Characteristics on the left have speed $1 > \frac12$ and catch up with the shock; on the right they have speed $0 < \frac12$ and the shock overtakes them.

Solution

For $t < 2$: the fan from $0$ is the Riemann solution for $0 \to 1$; the jump from $1$ to $0$ at $x = 1$ is a shock at speed $\frac12$. They meet when $t = 1 + \frac t2$, i.e. $t = 2$, at $x = 2$. Afterwards the left state at the shock is $\frac st$, so $\dot s = \frac12\cdot\frac st$; with $s(2) = 2$, $s = \sqrt{2t}$. Area: $\int_0^{\sqrt{2t}}\frac xt\,dx = \frac{2t}{2t} = 1$. Before $t = 2$ it is $\frac t2 + (1 + \frac t2 - t) = 1$ too. The height $\sqrt{2/t}$ decays like $t^{-1/2}$: the positive half of an N-wave.

Solution

With $\varphi(y, 0) = \exp(-\frac{1}{2\varepsilon}\int_0^yg)$, $\varphi = (4\pi\varepsilon t)^{-1/2}\int\exp\big(-\frac{(x - y)^2}{4\varepsilon t} - \frac{1}{2\varepsilon}\int_0^yg\big)dy$, and the exponent is $-\frac{1}{2\varepsilon}\big[\frac{(x - y)^2}{2t} + \int_0^yg\big]$. Differentiating in $x$ brings down $-\frac{x - y}{2\varepsilon t}$, so $-2\varepsilon\varphi_x$ is the same integral with the extra factor $\frac{x - y}{t}$. The normalising constants cancel in the ratio.

Solution

$F(\rho) = 100\rho(1 - \rho/150)$ cars/h. (a) $\sigma = \frac{F(45) - F(150)}{45 - 150} = \frac{100\cdot45\cdot0.7}{-105} = -30$ km/h: the back of the queue moves backward at $30$ km/h. (b) In $60$ s it moves $0.5$ km, so the queue is $0.5$ km long (about $75$ cars). (c) $F(75) = 100\cdot75\cdot0.5 = 3750$ cars/h.

Solution

$\sigma = \frac{1 - (-1)}{1 - (-1)} = 1$. $F'(\pm1) = 3$. On the right, characteristics have speed $3 > \sigma = 1$, so they move away from the shock: the shock is not compressive on that side.

Solution

(a) $\eta(u)_t + q(u)_x = uu_t + u^2u_x = u(u_t + uu_x) = 0$. (b) Split $\int_a^b$ at the shock $s(t)$ and differentiate: $\frac{d}{dt}\int_a^s\eta + \frac{d}{dt}\int_s^b\eta = \dot s\eta_l - \dot s\eta_r + \int\eta_t$, and $\int_a^s\eta_t + \int_s^b\eta_t = -\int q_x = q(a) - q_l + q_r - q(b)$. (c) With $d = u_r - u_l$: $-\sigma(\eta_r - \eta_l) = -\frac{(u_l + u_r)}{2}\cdot\frac{d(u_l + u_r)}{2}$ and $q_r - q_l = \frac{d(u_r^2 + u_ru_l + u_l^2)}{3}$. Their sum is $\frac{d}{12}\big[4(u_r^2 + u_ru_l + u_l^2) - 3(u_l + u_r)^2\big] = \frac{d}{12}(u_r - u_l)^2 = \frac{d^3}{12}$.