In [1]:

```
%matplotlib inline
```

In [2]:

```
%config InlineBackend.figure_format = 'svg'
from utils.snapshot_widgets import interact
from ipywidgets import widgets
from ipywidgets import FloatSlider, fixed
from exact_solvers import burgers
from exact_solvers import burgers_demos
from IPython.display import HTML
```

In this chapter, we study a simple scalar nonlinear conservation law: Burgers' equation. Burgers' equation models momentum transport in a fluid of uniform density and pressure, and it is the simplest equation that captures some key features of gas dynamics or water waves.

To examine the Python code for this chapter, see:

Burgers' equation has been used extensively for developing both theory and numerical methods, and it will allow us to explore the Riemann problem for a nonlinear conservation law. Burgers' equation is a scalar conservation law with flux $f(q)=q^2/2$:
\begin{align}
q_t + \left(\frac{1}{2}q^2\right)_x = 0.
\label{burgers0}
\end{align}
The quasilinear form is obtained by applying the chain rule to the flux term:

\begin{align*}
q_t + qq_x = 0.
\end{align*}

This equation looks very similar to the advection equation, with the difference that the advection speed at each point is given by the solution $q$. Burgers' equation is often viewed as a simplified version of equations in fluid dynamics or water waves that also have nonlinear fluxes. Our study of the dynamics of equation (\ref{burgers0}) is a step toward understanding these more complex nonlinear systems, in Shallow_water and Euler. In Traffic_flow we consider another scalar nonlinear conservation law that has a similar structure to Burgers' equation and a simpler physical interpretation.

The *characteristic speed* for Burgers' equation is $f'(q) = q$. As long as the solution is smooth, the solution is constant along characteristic curves $X(t)$ satisfying $X'(t) = f'(q(X(t),t)$ since then
$$
\frac{d}{dt} q(X(t),t) = q_x(X(t),t)X'(t) + q_t(X(t),t) = 0.
$$
Because of this, the characteristics are straight lines.
However, since the charactistic speed depends on the solution, these lines are not parallel and characterstics may converge or spread out.

In the figure below we consider Burgers' equation with a Gaussian hump as the initial data. Since the characteristic speed in Burgers' equation is given by $q$ itself, the peak of the hump travels faster than the rest, and characteristics are converging at the front of the traveling wave (where $f'(q)$ decreases with $x$) while they are spreading out behind the peak (where $f'(q)$ increases with $x$). The dashed line shows the initial condition while the solid lines show the solution at later times.

In [3]:

```
%matplotlib inline
```

In [4]:

```
%config InlineBackend.figure_format = 'svg'
from utils.snapshot_widgets import interact
from ipywidgets import widgets
from ipywidgets import FloatSlider, fixed
from exact_solvers import burgers
from exact_solvers import burgers_demos
```

In [5]:

```
interact(burgers_demos.multivalued_solution,
t=FloatSlider(min=0.,max=9.,value=7.75), fig=fixed(0));
```

Notice that at first $q$ remains single-valued for every $x$. However, after some time the crest of the wave overtakes the leading edge. After this time, we obtain a triple-valued solution for certain values of $x$. The first time this overtaking happens is referred to as the *breaking time* -- a reference to waves breaking on the beach. It is also the point where the conservation law, in its differential form, breaks down and where the characteristics cross each other for the first time. When characteristics cross, a shock wave, or discontinuity, forms. Mathematically, replacing the triple-valued region with a discontinuity will avoid the problem of the solution being multivalued. Where should the shock be located?

If we replace part of the multivalued solution interval with a shock, some mass will be removed (area $A_1$ in the figure below) and some mass will be added (area $A_2$ in the figure below). In the live notebook, the figure below shows possible locations for the shock at a given time; in order to maintain conservation of the integral of $q$, the shock must be placed so that these two areas are equal.

In [6]:

```
interact(burgers_demos.shock_location,
xshock = FloatSlider(min=6,max=10,step=0.25,value=7.75),
fig=fixed(0));
```

This geometric reasoning provides a nice intuition for the shock location, but is cumbersome in practice. To determine the location of the shock we can use the *Rankine-Hugoniot jump condition,* which we will derive in Traffic_flow, and which requires that the jump in flux across a propagating shock must be related to the jump in $q$ by
$$
s(q_r-q_l) = f(q_r) - f(q_l),
$$
at each instant in time, where $s$ is the shock speed at this time.

Since the flux for Burgers' equation is $f(q) = q^2/2$, this gives
\begin{align*}
s(q_r-q_l) & = \frac{1}{2} (q_r^2 - q_l^2) \\
\Rightarrow \ \ \ \ s &= \frac{1}{2}(q_l + q_r).
\end{align*}

In general (as in the image above) the states $q_l$ and $q_r$ just to the left and right of the shock will not be constant and the speed of a shock will change in time.

For Burgers' equation, (or any scalar hyperbolic conservation law with a *convex* flux function, such as the traffic flow problem considered in Traffic_flow), the solution of the Riemann problem consists of a single wave, which may be either a shock or a rarefaction, separating regions of constant value $q_l$ and $q_r$. Here convexity of the flux $f(q)$ means that $f'(q)$ is either monotonically increasing or monotonically decreasing as $q$ varies. If $f(q)$ is not convex then the solution can be more complicated; see the section on Convexity below and in Nonconvex_scalar.

The analysis above already yields the solution to the Riemann problem in the case that the resulting wave is a shock. The entropy condition, as we will explain below, indicates that a shock will form when $f'(q_l) > f'(q_r)$; i.e. when $q_l>q_r$. In this case, the solution is
\begin{align*}
q(x,t) =
\begin{cases}
q_l \quad \text{if} \ \ x/t<s \\
q_r \quad \text{if} \ \ x/t>s.
\end{cases}
\end{align*}

Below, we plot the solution of Burgers' equation for an initial condition that leads to a shock (i.e. with $q_l>q_r$).

In [7]:

```
interact(burgers_demos.shock(),
t=widgets.FloatSlider(min=0,max=1.0,value=0.5),
which_char=widgets.Checkbox(value=True,
description='Show characteristics'));
```

In the previous figure with the hump as the initial condition, we observed that a shock formed on the right side of the hump. However, on the left side, the characteristics spread out and will never cross. This part of the solution is called a rarefaction wave. This is the kind of behavior we will observe in the solution of the Riemann problem when $q_l<q_r$.

In the next figure, we consider such a Riemann problem. Although the initial data is discontinuous at $x=0$, we can think of smearing it out slightly to a continuous function so that all values between $q_l$ and $q_r$ are taken over a very narrow region around $x=0$. As time evolves, each value of $q$ propagates according to the advection speed $q$ given by the quasi-linear equation. Therefore, after time $t$, each value $q$ must propagate a distance $x=qt$, so the solution along the rarefaction is then $q = x/t$. As $q_l<q_r$, the smallest and largest displacements are given by $q_l t$ and $q_r t$, respectively.

In [8]:

```
interact(burgers_demos.rarefaction_figure,
t = FloatSlider(min=0,max=9,value=5));
```

The full rarefaction solution for the Burgers' Riemann problem is then simply given by \begin{align*} q(x,t) = \begin{cases} q_l, \quad \text{for} \ \ x<q_l t \\ x/t, \quad \text{for} \ \ q_l t \le x \le q_r t \\ q_r, \quad \text{for} \ \ x>q_r t. \end{cases} \end{align*}

As we will see in Traffic_flow, the rarefaction solution is always a self-similar solution. This means that it can be expressed as a function of the ratio between position and time $q(x,t) = \tilde{q}(x/t)$, so it remains the same when rescaling both $x$ and $t$ by the same factor. In Burgers' equation, the form of the rarefaction is particularly simple since the advection speed is simply $q$ and so the rarefaction wave is linear in $x$ with slope $1/t$ at time $t$.

In the figure below, we plot a solution of the Riemann problem with $q_l<q_r$, exhibiting a rarefaction.

In [9]:

```
interact(burgers_demos.rarefaction(),
t=widgets.FloatSlider(min=0,max=1.0,value=0.5),
which_char=widgets.Checkbox(value=True,
description='Show characteristics'));
```

As we mentioned before, the differential form of the equation breaks down in the presence of shocks/discontinuities. However, the integral form of the conservation law remains valid. Let's integrate the general conservation law $q_t+f(q)_x=0$ from $x=x_1$ to $x=x_2$ and $t=t_1$ to $t=t_2$: \begin{align*} \int_{t_1}^{t_2}\int_{x_1}^{x_2} [q_t+f(q)_x] dx dt = 0. \end{align*} This integral can be rewritten in terms of an indicator function $\phi(x,t)$ in $[x_1,x_2]\times[t_1,t_2]$ (defined to be 1 in this region, zero elsewhere): \begin{align} \int_{0}^{\infty}\int_{-\infty}^{\infty} [q_t+f(q)_x]\phi(x,t) dx dt = 0. \label{eq:Burgersintclaw2} \end{align}

We can further replace $\phi(x,t)$ by a smooth function with compact support on some region of the $x-t$ plane (zero outside a closed and bounded region). Assuming $t=0$ is in the support of $\phi(x,t)$, integration by parts yields
\begin{align}
\int_{0}^{\infty}\int_{-\infty}^{\infty} [q\phi_t+f(q)\phi_x] dx dt = -\int_{-\infty}^{\infty}q(x,0)\phi(x,0)dx,
\label{eq:Burgersintclaw3}
\end{align}
where now the derivatives are on $\phi(x,t)$ and not on $q$ or $f(q)$, so the equation still makes sense for discontinuous $q$. Note that we only obtain one boundary term along $t=0$ since $\phi(x,t)$ vanishes at infinity. The function $q(x,t)$ is called a *weak solution* of the conservation law if (\ref{eq:Burgersintclaw3}) holds for all functions $\phi(x,t)$ that are continuously differentiable and have compact support (bump functions). However, the function $\phi(x,t)$ chosen in (\ref{eq:Burgersintclaw2}) does not satisfy these conditions since it is not smooth. Nonetheless, we can approximate this function arbitrarily well by a smooth function. Note that any weak solution is a solution of the integral conservation law and vice versa.

Weak solutions are thus allowed to include discontinuities. The shock solution presented above, for instance, is a weak solution of Burgers' equation. After characteristics cross, there is no strong solution of the conservation law, and we must resort to weak solutions.

A potential problem with the notion of weak solutions is that in some cases the weak solution is not unique. For example, consider the Riemann problem for Burgers' equation with $q_l < q_r$. As mentioned already we expect a rarefaction to be the correct solution. This is because if we smeared out the initial data just slightly then there would be smoothly varying characteristics emerging from each point $x$ and these characteristics would spread out. However, for exactly discontinuous initial data, there exists another weak solution in which the initial discontinuity propagates as a shock wave with the Rankine-Hugoniot speed $(q_l + q_r)/2$: \begin{align*} q(x,t) = \begin{cases} q_l \quad \text{if} \ \ x/t<s \\ q_r \quad \text{if} \ \ x/t>s. \end{cases} \end{align*} This unphysical solution is also referred to as an expansion shock, and it is also a weak solution to Burgers' equation. In the next figure, we plot this weak solution.

In [10]:

```
interact(burgers_demos.unphysical(),
t=widgets.FloatSlider(min=0,max=1.0,value=0.5),
which_char=widgets.Checkbox(value=True,
description='Show characteristics'));
```

Note the behavior of the characteristics with respect to the shock. The fact that the characteristics are spreading away from the shock rather than converging on it indicates that the correct solution should instead be a rarefaction. In order to be able to specify which of the weak solutions is physically correct, we need to derive a mathematical condition from our physical intuition gained from observing the behavior of the characteristics. This condition is referred to as the entropy condition.

A given initial value problem for a hyperbolic PDE may have many weak solutions.
A condition that selects a unique physically correct solution out of these weak solutions is called an *admissibility condition* or more often an *entropy condition*. This name comes from gas dynamics, where the physical entropy must increase in the gas passing through a shock, according to the second law of thermodynamics. A discontinuity that violates this is non-physical. A mathematical "entropy function" with a similar property can often be defined for other conservation laws. More discussion of this and several other formulations of admissibility conditions, with more detailed explanations, are available in the literature, for example in many of the books cited in the Preface.

In the context of the Riemann problem, the entropy condition allows us to determine if the physical solution should involve a shock or a rarefaction. In the case of scalar conservation laws, the entropy condition is quite simple and can be formulated in terms of the flux function of the conservation law as follows: the solution of a scalar Riemann problem will consist of a shock only if
$$
f'(q_l) > f'(q_r).
$$
In other words, the solution is a shock only if nearby characteristics from the left and right approach the shock as time progresses. This is often called the *Lax Entropy Condition*. If nearby characteristics are spreading out, as in the last example, the correct solution is instead a rarefaction. In the case of Burgers' equation, the flux function is $f(q)=q^2/2$, so the correct solution is a shock only if
$$
q_l > q_r,
$$
which can be clearly observed in the interactive solution below (in the notebook). In later chapters, we will see how this condition generalizes to systems of conservation laws.

The flux $f(q) = q^2/2$ for Burgers' equation is a convex function, since $f''(q) = 1 > 0$ for all values of $q$. This means that as $q$ varies between $q_l$ and $q_r$ the the characteristic speed $f'(q) = q$ is either monotonically increasing (if $q_l < q_r$) or monotonically decreasing (if $q_l > q_r$). Hence for any Riemann problem the characteristic speeds for intermediate states are either purely converging or purely diverging, and the Riemann solution is always either a single rarefaction wave or a single shock wave.

The flux $f(\rho) = \rho(1-\rho)$ of the simple traffic flow model that we will consider in Traffic_flow also has the property that $f'(\rho)$ varies monotonically with $\rho$ (in the context of hyperbolic PDEs, a flux function is referred to as *convex* if either $f''(q)\ge0$ for all $q$ or $f''(q)\le 0$ for all $q$). Since $f''$ is always negative for the traffic flow model, the correct Riemann solution is a shock if $\rho_l < \rho_r$, or a rarefaction wave if $\rho_l > \rho_r$.

The solution to the Riemann problem can be much more complicated if $f'(q)$ is not monotonically varying between $q_l$ and $q_r$, i.e. if $f''(q)$ changes sign. In this case the Riemann solution can consist of multiple shock and rarefaction waves. The nonconvex case is explored further in Nonconvex_scalar.

In the live notebook, the figure below is an interactive solution of the Riemann problem for Burgers' equation. The values of the initial conditions and the time can be modified to observe their effect on the characteristic structure and the solution.

In [11]:

```
burgers.riemann_solution_interact()
```

We can now explore a few more examples that are representative of phenomena we will observe in more complicated systems, with more interesting initial data that the single jump discontinuity of a Riemann problem. For the animations shown below, the solution is computed numerically using PyClaw.

The animation below in the notebook, or on this webpage, shows the solution of Burgers' equation for an initial Gaussian hump, as discussed at the beginning of this chapter. A shock forms on the leading edge and the trailing edge spreads out as a rarefaction.

In [12]:

```
anim = burgers_demos.bump_animation(numframes = 50)
HTML(anim)
```

Out[12]: