Classification and first-order equations

Types of PDE, well-posedness and characteristics.

In this chapter 6

Tick a unit to add it to your PDF basket.

Classification and characteristics

Classification and characteristics

Definition

Classification of second-order equations

Adapted

Add to PDF

Defines elliptic, parabolic and hyperbolic second-order equations by the sign of the discriminant (and, in many variables, by the eigenvalue signs of the coefficient matrix), identifies the type with the number of real characteristic directions, reduces the two-variable constant-coefficient equation to Laplace, wave or heat canonical form, and states what well-posedness in Hadamard's sense asks for.

Prerequisites: Partial derivatives and the chain rule, Quadratic forms and the Rayleigh quotient, Spectral theorem for Hermitian matrices

Definition(The equation, and how nonlinear it is)AdaptedAdapted from Evans (2010), §1.1, Definitions (i)–(iv), p. 2; Arfken, Weber & Harris (2013), §9.3, Eq. (9.17), p. 410

Throughout, $u$ is the unknown and subscripts denote partial derivatives. A partial differential equation of order $k$ is linear if it is a sum of the derivatives of $u$ of order at most $k$ with coefficients depending on the independent variables alone, equal to a given function of those variables; semilinear if its order- $k$ terms have that form and the remaining terms may depend on $u$ and its lower derivatives; quasilinear if the coefficients of the order- $k$ derivatives may depend on $u$ and its derivatives of order up to $k-1$ ; and fully nonlinear if it depends nonlinearly on the derivatives of order $k$ themselves.

Classification concerns the order- $k$ terms only, so the whole hierarchy up to quasilinear is covered by one computation. For second order in two variables we write

a(x,y)\,u_{xx} + 2b(x,y)\,u_{xy} + c(x,y)\,u_{yy} = F\!\left(x, y, u, u_x, u_y\right), \tag{1}

with $a, b, c$ real and not all zero at the point under discussion. The factor $2$ on the mixed term is a convention; it makes the second-order part the quadratic form of the symmetric matrix $\begin{pmatrix} a & b \\ b & c\end{pmatrix}$ . The discriminant of (1) at a point is

\delta = b^2 - ac .

Riley, Hobson and Bence instead write $A u_{xx} + B u_{xy} + C u_{yy}$ , so that their discriminant is $B^2 - 4AC$ ; the two agree because $B = 2b$ , and all statements below are unchanged if $\delta$ is replaced by $B^2 - 4AC$ .

Definition(Elliptic, parabolic, hyperbolic)AdaptedAdapted from Stakgold & Holst (2011), §8.1, pp. 465–466; Strauss (2008), §1.6, Thm 1 and Example 1, pp. 28–29

Equation (1) is called, at the point $P$ ,

elliptic if $\delta < 0$ at $P$ ,
parabolic if $\delta = 0$ at $P$ ,
hyperbolic if $\delta > 0$ at $P$ .

It is called elliptic (parabolic, hyperbolic) on a region if it is so at each of its points. The names come from the conic $a\xi^2 + 2b\xi\eta + c\eta^2 = \text{const}$ , which is an ellipse, a parabola or a hyperbola in the same three cases.

The prototypes, with $(x,y)$ read as $(x,t)$ where a time appears, are Laplace's equation $u_{xx} + u_{yy} = 0$ ( $\delta = -1$ , elliptic), the wave equation $u_{tt} - c^2 u_{xx} = 0$ ( $\delta = c^2 > 0$ , hyperbolic) and the diffusion equation $u_t - k u_{xx} = 0$ ( $\delta = 0$ , parabolic, since there is no $u_{tt}$ term and no mixed term).

Definition(The type counts the real characteristic directions)AdaptedAdapted from Stakgold & Holst (2011), §8.1, Eqs. (8.1.14)–(8.1.15) and the three cases, pp. 465–466; Hassani (2013), Example 21.1.6, Eq. (21.7), pp. 641–642

Let $\Gamma$ be a smooth curve in the $x$ – $y$ plane, let $P \in \Gamma$ , and suppose that $u$ and both of its first derivatives are prescribed along $\Gamma$ (Cauchy data). Call $\Gamma$ characteristic at $P$ if these data together with equation (1) fail to determine all three second derivatives $u_{xx}, u_{xy}, u_{yy}$ at $P$ . Then $\Gamma$ is characteristic at $P$ if and only if its direction $(dx, dy)$ there satisfies

a\,dy^2 - 2b\,dx\,dy + c\,dx^2 = 0 . \tag{2}

Consequently, if $a(P) \ne 0$ , the characteristic directions at $P$ are the solutions of

\frac{dy}{dx} = \frac{b \pm \sqrt{b^2 - ac}}{a} , \tag{3}

and there are two of them when $\delta > 0$ , one when $\delta = 0$ , and none (in the reals) when $\delta < 0$ . A curve characteristic at each of its points is a characteristic curve; where the equation is hyperbolic, integrating the two branches of (3) gives two families of characteristic curves, and where it is parabolic, one family.

Idea of the proofOriginal

Cauchy data on $\Gamma$ give $u_x$ and $u_y$ along $\Gamma$ , so differentiating those two functions along $\Gamma$ gives two linear relations among the three second derivatives at $P$ . The equation itself gives a third. Three linear equations in three unknowns are solvable exactly when their determinant does not vanish, and that determinant is the quadratic form in $(dx, dy)$ appearing in (2). So the characteristic directions are the null directions of that form, and how many there are is decided by the sign of its discriminant.

ProofAdaptedAdapted from Hassani (2013), Example 21.1.6, pp. 641–642; Stakgold & Holst (2011), Eq. (8.1.14), p. 466

Write $p = u_x$ and $q = u_y$ , both known on $\Gamma$ by hypothesis, and let $(dx, dy)$ be a tangent direction to $\Gamma$ at $P$ .

Step 1 (two relations from the data). Moving along $\Gamma$ by the infinitesimal step $(dx, dy)$ changes $p$ and $q$ by amounts the data already know. Since $p$ and $q$ are differentiable, the chain rule gives

dp = u_{xx}\,dx + u_{xy}\,dy, \qquad dq = u_{xy}\,dx + u_{yy}\,dy ,

with $dp, dq$ computable from the prescribed $p, q$ along $\Gamma$ .

Step 2 (one relation from the equation). Evaluating (1) at $P$ and using that $u, u_x, u_y$ are known there, the right-hand side $F$ is a known number $-\!f$ , so

a\,u_{xx} + 2b\,u_{xy} + c\,u_{yy} = -f .

Step 3 (solve or not). Steps 1 and 2 are three linear equations for $(u_{xx}, u_{xy}, u_{yy})$ with coefficient matrix

M = \begin{pmatrix} dx & dy & 0 \\ 0 & dx & dy \\ a & 2b & c \end{pmatrix}, \qquad \det M = c\,dx^2 - 2b\,dx\,dy + a\,dy^2 .

(The determinant is computed by expanding along the first row.) The system has a unique solution precisely when $\det M \ne 0$ , so by the definition of "characteristic at $P$ ", $\Gamma$ is characteristic at $P$ exactly when $\det M = 0$ , which is (2).

Step 4 (count the directions). A direction is fixed up to scale, so put $m = dy/dx$ when $dx \ne 0$ . Dividing (2) by $dx^2$ gives $a m^2 - 2bm + c = 0$ . If $a \ne 0$ this is a genuine quadratic with discriminant $4(b^2 - ac) = 4\delta$ , whence the formula (3) and the count two, one or none according to $\delta > 0$ , $\delta = 0$ , $\delta < 0$ . If $a = 0$ then $\delta = -ac + b^2 = b^2 \ge 0$ , so the elliptic case cannot occur, and (2) becomes $dy\,(-2b\,dx + c\,dy) = 0$ : the direction $dy = 0$ is characteristic, and it is the only one when $b = 0$ (then $\delta = 0$ ), while a second direction $dy/dx = 2b/c$ appears when $b \ne 0$ and $c \ne 0$ (then $\delta = b^2 > 0$ ). The count is as stated in both cases. $\square$

RemarkOriginal

Where the hypotheses are used, and what the count means.

Cauchy data, not just $u$ . Step 1 needs $u_x$ and $u_y$ along $\Gamma$ . If only $u$ is given, there is one relation instead of two and the argument says nothing; the Dirichlet problem is not a Cauchy problem.
Three unknowns, three equations. The determinant condition is what makes the classification a statement about the second-order coefficients only. The lower-order terms enter only through the known number $f$ in Step 3, which is why $F$ may depend on $u, u_x, u_y$ in any way at all, linearly or not.
Quasilinear equations. If $a, b, c$ depend on $u, u_x, u_y$ as well as on $x, y$ , the argument still works at a point once a solution is fixed there, because those arguments are known from the data. The type can then differ from solution to solution.
Why "none" is not the same as "harmless". The absence of real characteristics means no direction is singled out for propagation, which is why elliptic equations have no wave-like behaviour. It does not mean the Cauchy problem for an elliptic equation is solvable; the example below shows it is not.

Definition(Reduction to canonical form)AdaptedAdapted from Strauss (2008), §1.6, Thm 1, pp. 28–29

Write the two-variable linear equation with six real constant coefficients as

a_{11}u_{xx} + 2a_{12}u_{xy} + a_{22}u_{yy} + a_1 u_x + a_2 u_y + a_0 u = 0 , \tag{4}

whose discriminant is $\delta = a_{12}^2 - a_{11}a_{22}$ , and exclude the degenerate case in which $a_{11}, a_{12}, a_{22}$ all vanish. Then some invertible linear change of the independent variables, $(x,y) \mapsto (\xi,\eta)$ , turns the second-order part of (4) into exactly

u_{\xi\xi} + u_{\eta\eta} \quad (\delta < 0), \qquad u_{\xi\xi} - u_{\eta\eta} \quad (\delta > 0), \qquad u_{\xi\xi} \quad (\delta = 0) ,

up to an overall nonzero factor, while the terms of order one and zero remain terms of order one and zero. Strauss records this as the reduction of (4) to $u_{xx} + u_{yy} + \cdots = 0$ , $u_{xx} - u_{yy} + \cdots = 0$ or $u_{xx} + \cdots = 0$ , with the dots standing for the lower-order terms, and calls these the elliptic, hyperbolic and parabolic cases respectively. So each of Laplace's equation, the wave equation and the one-dimensional diffusion equation is not merely an instance of its class but a normal form for the whole class.

Idea of the proofOriginal

The second-order part is a quadratic form in the two differentiation operators. Completing the square turns it into $(\partial_x + a_{12}\partial_y)^2$ plus a multiple of $\partial_y^2$ , and the sign of that multiple is the sign of $-\delta$ . A linear change of variables can be chosen so that the squared operator becomes $\partial_\xi$ and the leftover $\partial_y$ becomes $\partial_\eta$ , which is exactly the three canonical forms. The only work beyond the square is arranging $a_{11} \ne 0$ first, which costs at most a swap of the variables or one preliminary rotation.

ProofAdaptedAdapted from Strauss (2008), §1.6, proof of Thm 1, Eqs. (2)–(4), p. 29

Linear changes of variable act on the second-order coefficient matrix $A = \begin{pmatrix} a_{11} & a_{12} \\ a_{12} & a_{22}\end{pmatrix}$ by $A \mapsto BAB^{\mathsf T}$ , so they multiply $\det A = a_{11}a_{22} - a_{12}^2 = -\delta$ by $(\det B)^2 > 0$ and cannot change the sign of $\delta$ . Hence the three cases are mutually exclusive and preserved, and we may normalise freely.

Step 1 (make $a_{11} = 1$ ). Suppose first that some diagonal coefficient is nonzero. If $a_{11} = 0$ and $a_{22} \ne 0$ , interchange $x$ and $y$ ; this is linear and leaves $a_{12}$ and $\delta$ unchanged. So assume $a_{11} \ne 0$ and divide the whole equation by $a_{11}$ : the new coefficients are $1$ , $a_{12}/a_{11}$ , $a_{22}/a_{11}$ , and the new discriminant is $\delta/a_{11}^2$ , of the same sign. Rename them $1, a_{12}, a_{22}$ , so now $\delta = a_{12}^2 - a_{22}$ .

Step 2 (complete the square). With $\partial_x, \partial_y$ written as operators,

\partial_x^2 + 2a_{12}\partial_x\partial_y + a_{22}\partial_y^2 = \bigl(\partial_x + a_{12}\partial_y\bigr)^2 + \bigl(a_{22} - a_{12}^2\bigr)\partial_y^2 = \bigl(\partial_x + a_{12}\partial_y\bigr)^2 - \delta\,\partial_y^2 ,

as expanding the square confirms. The constant-coefficient operators commute, so no ordering question arises.

Step 3 (the elliptic case $\delta < 0$ ). Put $\beta = \sqrt{-\delta} = (a_{22} - a_{12}^2)^{1/2} > 0$ and change variables by

x = \xi, \qquad y = a_{12}\,\xi + \beta\,\eta ,

which is invertible since $\beta \ne 0$ . By the chain rule, $\partial_\xi = \partial_x + a_{12}\partial_y$ and $\partial_\eta = \beta\,\partial_y$ . Substituting into Step 2 turns the second-order part into $\partial_\xi^2 + \beta^2 \partial_y^2 = \partial_\xi^2 + \partial_\eta^2$ , Laplace's operator. The first-order and zeroth-order terms become first-order and zeroth-order terms in $\xi, \eta$ , because the transformation is linear.

Step 4 (the hyperbolic case $\delta > 0$ with $a_{11} \ne 0$ ). Now put $\beta = \sqrt{\delta} > 0$ and use the same change of variables $x = \xi$ , $y = a_{12}\xi + \beta\eta$ . Step 2 becomes $\partial_\xi^2 - \beta^2\partial_y^2 = \partial_\xi^2 - \partial_\eta^2$ , the wave operator.

Step 5 (the parabolic case $\delta = 0$ ). Here Step 2 leaves only $(\partial_x + a_{12}\partial_y)^2$ . Take $x = \xi$ , $y = a_{12}\xi + \eta$ , so that $\partial_\xi = \partial_x + a_{12}\partial_y$ and $\partial_\eta = \partial_y$ ; the second-order part is $\partial_\xi^2$ , and the remaining terms are of order $1$ or $0$ . This is case 3.

Step 6 (both diagonal coefficients zero). If $a_{11} = a_{22} = 0$ then $\delta = a_{12}^2$ , and the excluded case $a_{11} = a_{12} = a_{22} = 0$ aside, $a_{12} \ne 0$ , so $\delta > 0$ : the equation must be hyperbolic, and Step 1 is unavailable. Take $x = \xi + \eta$ , $y = \xi - \eta$ , so $\partial_\xi = \partial_x + \partial_y$ and $\partial_\eta = \partial_x - \partial_y$ , whence

\partial_\xi^2 - \partial_\eta^2 = (\partial_x + \partial_y)^2 - (\partial_x - \partial_y)^2 = 4\,\partial_x\partial_y .

So $2a_{12}u_{xy} = \tfrac{a_{12}}{2}\bigl(u_{\xi\xi} - u_{\eta\eta}\bigr)$ , and dividing by $a_{12}/2$ gives case 2. $\square$

RemarkOriginal

What the reduction does and does not give. Steps 3 to 5 use a single linear map, so for constant coefficients one change of variables works on the whole plane. If $a, b, c$ vary with $x$ and $y$ , the same algebra works at each point separately, so the type is still well defined pointwise, but the transformation then varies from point to point and one must integrate (3) to produce coordinates; the canonical form is achieved only locally, and only where the type does not change. Arfken reaches the same conclusion by the explicit substitution $\xi = c^{1/2}x - c^{-1/2}by$ , $\eta = c^{-1/2}y$ , which requires $c > 0$ , and obtains $L = (ac - b^2)\partial_\xi^2 + \partial_\eta^2$ ; the coefficient $ac - b^2 = -\delta$ again carries the sign that decides the type.

Definition(Many variables)AdaptedAdapted from Strauss (2008), §1.6, pp. 29–31

For $n$ variables $x_1, \dots, x_n$ and a linear equation

\sum_{i,j=1}^{n} a_{ij}\,u_{x_i x_j} + \sum_{i=1}^{n} a_i\,u_{x_i} + a_0 u = 0 \tag{5}

with real constants and $a_{ij} = a_{ji}$ , let $d_1, \dots, d_n$ be the eigenvalues of $A = (a_{ij})$ . A linear change of variables $\xi = Bx$ sends $A$ to $BAB^{\mathsf T}$ , and by the spectral theorem for real symmetric matrices there is a rotation $B$ making this diagonal; a further change of scale replaces each nonzero diagonal entry by $\pm 1$ . Equation (5) is called

elliptic if $d_1, \dots, d_n$ are either all positive or all negative — equivalently, if $A$ or $-A$ is positive definite;
hyperbolic if $0$ is not among them and the sign pattern is $n-1$ of one sign against a single one of the other;
ultrahyperbolic if $0$ is not among them and each sign occurs two or more times;
parabolic if $0$ occurs exactly once among them and one sign accounts for the remaining $n-1$ .

With $x_n = t$ , Laplace's and Poisson's equations are elliptic, the wave and Klein-Gordon equations are hyperbolic, and the diffusion equation is parabolic. Ultrahyperbolic equations occur rarely in physics and are not pursued here. For variable coefficients the same test is applied at each point, and the type may differ from region to region.

RemarkAdaptedAdapted from Hassani (2013), p. 641; Evans (2010), §7.1.1, Definition and Remark, p. 350

Two things worth keeping straight.

The classification is a statement about a signature. Diagonalising $A$ by $A \mapsto BAB^{\mathsf T}$ does not preserve the eigenvalues, only how many are positive, negative and zero; that is Sylvester's law of inertia (see quadratic forms). The four cases above are phrased only in terms of those counts, so they are invariant under every invertible linear change of variables, which is what makes the type a property of the equation.

Sources differ on "parabolic". Hassani, writing the equation with $(A_{jk})$ already diagonal as $\sum_j a_j u_{x_jx_j} + F = 0$ , calls it elliptic at $x_0$ if all $a_j(x_0)$ are nonzero with the same sign, hyperbolic if all are nonzero with exactly one of opposite sign, ultrahyperbolic if all are nonzero but with at least two of each sign, and parabolic if at least one $a_j(x_0)$ vanishes. That is weaker than the definition used here, which is Strauss's: it also calls parabolic, for instance, an equation whose coefficient matrix has one zero eigenvalue and both signs among the rest. This unit follows Strauss throughout; readers comparing with Hassani should check which convention a given statement needs. The two agree on all of the standard equations of physics listed above.

A third formulation, standard in the modern theory, replaces signs by a uniform bound: Evans calls $\partial_t + L$ with $Lu = -\sum a_{ij}(x,t)u_{x_ix_j} + \sum b_i u_{x_i} + cu$ uniformly parabolic if there is a constant $\theta > 0$ with $\sum_{i,j} a_{ij}(x,t)\xi_i\xi_j \ge \theta|\xi|^2$ for all $(x,t)$ in the region and all $\xi \in \mathbb{R}^n$ , so that $L$ is uniformly elliptic in $x$ at each fixed time. This is the pointwise condition "all eigenvalues positive" made quantitative and uniform, which is what the estimates of the modern theory need.

Definition(Well-posedness in Hadamard's sense)AdaptedAdapted from Strauss (2008), §1.5, pp. 25–26

A problem consisting of a partial differential equation in a domain together with initial and/or boundary conditions is well posed if

existence: there is at least one solution satisfying all the conditions;
uniqueness: there is at most one;
stability: the solution depends on the data in a stable way, in that a small change in the data produces only a small change in the solution.

Condition 3 presupposes a notion of nearness, that is, a norm or metric on the data and on the solutions, and a problem can be well posed for one choice and ill posed for another; which choice is intended must always be said. Too few auxiliary conditions and the problem is underdetermined, with many solutions; too many and it is overdetermined, with none. Stability is not a mathematical luxury: measured data carry errors, and a problem that magnifies them without bound cannot be solved from measurements at all.

Example(Cauchy data for Laplace's equation: no stability)AdaptedAdapted from Strauss (2008), §1.5, Eq. (2), p. 27

Take Laplace's equation $u_{xx} + u_{yy} = 0$ in the upper half-plane $D = \{-\infty < x < \infty,\ y > 0\}$ , and prescribe both $u$ and $u_y$ on the boundary $y = 0$ . Laplace's equation is elliptic, so by the statement above it has no real characteristics anywhere, and the obstruction of a characteristic initial curve is absent. Nevertheless the problem is not well posed. For $n = 1, 2, 3, \dots$ put

u_n(x,y) = \frac{1}{n}\,e^{-\sqrt{n}}\,\sin(nx)\,\sinh(ny) .

Each $u_n$ is a solution. Differentiating, $\partial_x^2 u_n = -n^2 u_n$ and $\partial_y^2 u_n = +n^2 u_n$ , so $u_{n,xx} + u_{n,yy} = 0$ in $D$ .
The data go to zero. On $y = 0$ , $\sinh 0 = 0$ gives $u_n(x,0) = 0$ , and $\partial_y u_n(x,y) = e^{-\sqrt n}\sin(nx)\cosh(ny)$ gives $\partial_y u_n(x,0) = e^{-\sqrt n}\sin(nx)$ . Hence $\sup_x |u_n(x,0)| = 0$ and $\sup_x|\partial_y u_n(x,0)| = e^{-\sqrt n} \to 0$ as $n \to \infty$ , and the same is true of every derivative of the data in $x$ , since differentiating $k$ times multiplies the bound by $n^k$ and $n^k e^{-\sqrt n} \to 0$ .
The solutions do not. Fix any $y_0 > 0$ . Since $\sinh(ny_0) \ge \tfrac12 e^{ny_0 - \log 2}$ for $n y_0 \ge 1$ , at the points $x = \pi/(2n)$ we get
$\bigl|u_n\bigl(\tfrac{\pi}{2n}, y_0\bigr)\bigr| = \frac{e^{-\sqrt n}\sinh(ny_0)}{n} \ \ge\ \frac{e^{\,n y_0 - \sqrt n}}{4n} \longrightarrow \infty .$
(The inequality for $\sinh$ holds because $\sinh z \ge e^z/4$ as soon as $e^{2z} \ge 2$ , hence for $z \ge 1$ .) For $y_0 = 1$ and $n = 25$ the bound reads $e^{20}/100 \approx 4.9 \times 10^6$ , and the exact value is $e^{-5}\sinh 25/25 \approx 9.70 \times 10^6$ , while the data are of size $e^{-5} \approx 6.7\times 10^{-3}$ .
Conclusion. Data arbitrarily close to zero, whose own solution is $u \equiv 0$ , produce solutions that are arbitrarily large at every interior height. Stability fails in the supremum norm, so the Cauchy problem for Laplace's equation is ill posed. Compare this with the hyperbolic case, where Cauchy data on a non-characteristic curve do give a unique, stably dependent solution; and note that for elliptic equations it is the boundary value problem, with $u$ or $\partial u/\partial n$ given on a closed boundary, that is well posed instead.

ExampleOriginal

(a) Classify and reduce $3u_{xx} + 10u_{xy} + 3u_{yy} = 0$ .

Read off the coefficients in the convention of (1). $a = 3$ , $2b = 10$ so $b = 5$ , and $c = 3$ . Then $\delta = b^2 - ac = 25 - 9 = 16 > 0$ : hyperbolic in the whole plane.
Characteristics. By (3), $dy/dx = (5 \pm 4)/3$ , that is $3$ or $\tfrac13$ . Integrating, the two families are $y - 3x = \text{const}$ and $3y - x = \text{const}$ .
Canonical form. Divide by $a_{11} = 3$ to get $u_{xx} + \tfrac{10}{3}u_{xy} + u_{yy} = 0$ , so $a_{12} = \tfrac53$ , $a_{22} = 1$ and $\delta' = \tfrac{25}{9} - 1 = \tfrac{16}{9}$ , with $\beta = \sqrt{\delta'} = \tfrac43$ . Step 4 of the proof prescribes $x = \xi$ , $y = \tfrac53\xi + \tfrac43\eta$ , that is $\xi = x$ , $\eta = \tfrac{3y - 5x}{4}$ . Substituting gives $3u_{xx} + 10u_{xy} + 3u_{yy} = 3\bigl(u_{\xi\xi} - u_{\eta\eta}\bigr)$ , so the equation becomes the wave equation $u_{\xi\xi} - u_{\eta\eta} = 0$ .
Check against the factorisation. The operator factors as $3\partial_x^2 + 10\partial_x\partial_y + 3\partial_y^2 = (3\partial_x + \partial_y)(\partial_x + 3\partial_y)$ , so $u = f(x - 3y) + g(3x - y)$ solves it for all twice-differentiable $f, g$ ; direct substitution confirms this. The level sets of $x - 3y$ have slope $\tfrac13$ and those of $3x - y$ have slope $3$ : the arbitrary functions travel along the characteristics found in step 2, as they must.

(b) An equation of mixed type. Take the equation of Tricomi type

u_{xx} + y\,u_{yy} = 0 ,

for which $a = 1$ , $b = 0$ and $c = y$ , so that

\delta = b^2 - ac = -y .

Type. The equation is hyperbolic in the open half-plane $y < 0$ , parabolic on the line $y = 0$ , and elliptic in $y > 0$ — for every $x$ . One equation thus carries both behaviours, separated by a line on which it degenerates.
Characteristics where they exist. Equation (2) reads $dy^2 + y\,dx^2 = 0$ , so for $y < 0$ the two real slopes are $dy/dx = \pm\sqrt{-y}$ , which is also what (3) gives since $a = 1$ . Separating variables, $dy/\sqrt{-y} = \pm\,dx$ integrates to $-2\sqrt{-y} = \pm x + \text{const}$ , so the two families of characteristic curves are
$x \pm 2\sqrt{-y} = \text{const} \qquad (y < 0) ,$
two families of parabolas $y = -\tfrac14(x - \text{const})^2$ tangent to the line $y = 0$ from below. Differentiating this parabola gives $dy/dx = (\text{const} - x)/2 = \pm\sqrt{-y}$ , confirming the slopes.
What happens at the degeneracy. The two slopes $\pm\sqrt{-y}$ coalesce as $y \to 0^-$ , so the two families merge into the single family $dy = 0$ predicted on the parabolic line, and every characteristic is tangent to that line. A problem posed across $y = 0$ needs data appropriate to a hyperbolic equation below it and to an elliptic one above it; no single prescription serves both, which is what makes equations of mixed type hard.

ApplicationOriginal

Steady compressible flow, and why the type changes at Mach 1. The linearised potential equation for steady two-dimensional flow past a thin body at free-stream Mach number $M_\infty$ is

(1 - M_\infty^2)\,\phi_{xx} + \phi_{yy} = 0 ,

with $\phi$ the perturbation potential and $x$ along the flow. Here $a = 1 - M_\infty^2$ , $b = 0$ , $c = 1$ , so $\delta = M_\infty^2 - 1$ . Subsonic flow ( $M_\infty < 1$ ) is elliptic: there are no characteristics, every point of the flow feels every part of the body, and the computation needs data on a closed boundary around the whole field. Supersonic flow ( $M_\infty > 1$ ) is hyperbolic, with characteristic slopes $dy/dx = \mp 1/\sqrt{M_\infty^2 - 1}$ from (3); these are the Mach lines, inclined at the Mach angle $\mu = \arctan\bigl(1/\sqrt{M_\infty^2-1}\bigr)$ , and disturbances travel only along them, which is why a supersonic body influences only the wedge downstream of it. At $M_\infty = 1$ the equation is parabolic and both descriptions break down; transonic flow is modelled instead by equations of mixed type, whose coefficients change sign inside the domain, exactly as in part (b) of the example above.

The same classification governs practical choices elsewhere. It decides which data a solver must be given — Cauchy data on a non-characteristic surface for hyperbolic problems (seismic wave propagation, electromagnetics in the time domain), boundary data on a closed surface for elliptic ones (electrostatics, steady heat conduction, incompressible pressure fields) — and which are hopeless: running diffusion backwards in time, or recovering a source from potential measurements on a surface, are precisely the ill-posed problems whose instability the Laplace example above exhibits, and they are attacked in practice only by regularisation, which deliberately replaces the problem by a nearby well-posed one.

Method

First-order equations and characteristics

Adapted

Add to PDF

Turns a first-order quasilinear equation with data on a curve into a system of ordinary differential equations along characteristics, gives the recipe and the non-characteristic condition under which it determines a unique local solution, proves it, and works out the transport equation and three further examples.

Prerequisites: Existence and uniqueness of solutions, Implicit and inverse function theorems, Partial derivatives and the chain rule

Definition(Quasilinear first-order equations and their characteristics)AdaptedAdapted from Evans (2010), §3.2.2b, Eq. (20), p. 97, and §3.2.5b, Eq. (54), p. 107 (the quasilinear form); Stakgold & Holst (2011), §8.1, Eq. (8.1.8) with the linear and semilinear cases, p. 463, and Eqs. (8.1.9)–(8.1.10), p. 464 (the determinant condition)

A first-order equation in two independent variables is quasilinear if it can be written

a(x,y,u)\,u_x + b(x,y,u)\,u_y = c(x,y,u) , \tag{1}

so that the derivatives of $u$ appear linearly but their coefficients may involve $u$ . This generality is Evans's: his Eq. (20) is exactly the case in which $u$ may enter the coefficient of every derivative. Stakgold and Holst's first-order equation (8.1.8) is the narrower $a(x,y)u_x + b(x,y)u_y + c(x,y,u) = 0$ , in which $a$ and $b$ do not depend on $u$ at all; the two narrower cases named next are theirs. Equation (1) is semilinear if $a$ and $b$ do not depend on $u$ but $c$ does, and linear if in addition $c = -q(x,y)u + F(x,y)$ .

In $n$ variables Evans writes the quasilinear equation with its source term on the left,

\mathbf{b}(x,u)\cdot Du + c(x,u) = 0 , \qquad x \in \mathbb{R}^n ,

where $\mathbf{b}$ is the whole coefficient vector of the gradient — not the scalar $b$ of (1) — and consequently his $c$ is the negative of the right-hand side of (1).

The characteristic equations of (1) are the ordinary differential equations

\frac{dx}{ds} = a(x,y,z), \qquad \frac{dy}{ds} = b(x,y,z), \qquad \frac{dz}{ds} = c(x,y,z) , \tag{2}

for three unknown functions $x(s), y(s), z(s)$ . A solution of (2) is a characteristic strip — a curve in $(x,y,z)$ -space, carrying both a point of the plane and the value of $u$ there. Its shadow $s \mapsto (x(s), y(s))$ in the plane of the independent variables is a projected characteristic, or simply a characteristic curve of (1). (Stakgold and Holst use the word strip for a different object, a region of the plane, and put it in quotation marks; the statement below records that usage. Bolded, the term always means the solution of (2).) The whole point of (2) is that the third equation closes the system: no derivative of $u$ appears in it, so the three functions can be computed without knowing $u$ anywhere else.

A smooth curve $\Gamma$ in the $x$ – $y$ plane carrying the data $u = g$ is called characteristic at a point $P \in \Gamma$ if $u_x$ and $u_y$ cannot be computed at $P$ from the data on $\Gamma$ together with equation (1). Writing $(dx, dy)$ for a tangent direction of $\Gamma$ at $P$ , differentiating $u$ along $\Gamma$ gives $u_x\,dx + u_y\,dy = du$ , with $du$ known; equation (1) gives $a u_x + b u_y = c$ , with $c$ known once $u(P)$ is. These two linear equations for $(u_x, u_y)$ are solvable precisely when

\begin{vmatrix} dx & dy \\ a & b \end{vmatrix} = b\,dx - a\,dy \ne 0 , \tag{3}

so $\Gamma$ is characteristic at $P$ exactly when $b\,dx - a\,dy = 0$ there, that is, when the vector $(a,b)$ evaluated at $P$ is tangent to $\Gamma$ . Where $a \ne 0$ this says the characteristic direction has slope $dy/dx = b/a$ , which is (2) with $s$ eliminated. A curve on which (3) holds at every point is called non-characteristic.

Method(Solving a quasilinear first-order equation by characteristics)AdaptedAdapted from Evans (2010), §3.2.2b and §3.2.4, pp. 97–98 and 102–105; Arfken, Weber & Harris (2013), §9.2, Eqs. (9.8)–(9.9), pp. 406–407

Given equation (1) and data $u = g$ on a curve $\Gamma$ :

Parametrise the data. Write $\Gamma$ as $\bigl(f(\tau), h(\tau)\bigr)$ with $f, h \in C^1$ and $(f', h') \ne (0,0)$ , and let $g(\tau)$ be the prescribed value of $u$ at the parameter $\tau$ .
Check that $\Gamma$ is non-characteristic. Verify
$J(\tau) \;=\; a\bigl(f, h, g\bigr)\,h'(\tau) - b\bigl(f, h, g\bigr)\,f'(\tau) \;\ne\; 0 ,$
which is condition (3) with $(dx, dy) = (f', h')$ , up to sign. If $J$ vanishes at a point, stop: the method gives nothing there.
Solve the characteristic equations. Integrate (2) in $s$ for each fixed $\tau$ , with initial conditions
$x(0,\tau) = f(\tau), \qquad y(0,\tau) = h(\tau), \qquad z(0,\tau) = g(\tau) .$
For a linear or semilinear equation the first two equations decouple from the third and are solved first; for a genuinely quasilinear equation all three must be solved together.
Invert the change of variables. Solve $x = x(s,\tau)$ , $y = y(s,\tau)$ for $s = s(x,y)$ and $\tau = \tau(x,y)$ . Step 2 guarantees this is possible near $\Gamma$ .
Read off the solution. Set
$u(x,y) = z\bigl(s(x,y), \tau(x,y)\bigr) .$
Say where it holds. The formula is valid on the region swept out by the characteristics on which step 4 succeeded. Check the extent of that region: in the quasilinear case characteristics may cross, and the formula fails there even though every step above was legitimate.

A useful shortcut for linear equations with constant $a, b$ : changing to the variables $s' = ax + by$ and $t' = bx - ay$ converts (1) into $(a^2+b^2)\,\partial u/\partial s' + \dots = 0$ , an ordinary differential equation in $s'$ with $t'$ as a parameter, which is steps 3 to 5 carried out in one move.

Method(Local existence and uniqueness)AdaptedAdapted from Evans (2010), §3.2.1, p. 94 (the smoothness hypothesis), §3.2.4, Lemma 2 and Thm 2, pp. 102–105, and §3.2.5b, pp. 107–108; Stakgold & Holst (2011), §8.1, p. 464

Evans fixes a standing hypothesis at the start of his treatment of characteristics: the equation and the data are smooth. In the notation of (1) that means $a, b, c$ smooth functions of $(x,y,z)$ on an open set, $\Gamma$ a smooth curve $\tau \mapsto (f(\tau), h(\tau))$ , and $g$ smooth. Suppose that at $\tau = \tau_0$

a\,h'(\tau_0) - b\,f'(\tau_0) \ne 0, \qquad a, b \text{ evaluated at } \bigl(f(\tau_0), h(\tau_0), g(\tau_0)\bigr) .

This is Evans's non-characteristic condition $\mathbf{b}(x^0, z_0)\cdot\boldsymbol{\nu}(x^0) \ne 0$ , with $z_0 = g(\tau_0)$ and $\boldsymbol{\nu}$ his outward unit normal to $\partial U$ at the data point, written out for a parametrised plane curve: the tangent of $\Gamma$ there is $(f', h')$ , so a normal is $(h', -f')$ and $(a,b)\cdot(h',-f') = a h' - b f'$ . Then there is a neighbourhood $V$ of the point $P_0 = (f(\tau_0), h(\tau_0))$ in which the problem

a(x,y,u)u_x + b(x,y,u)u_y = c(x,y,u), \qquad u = g \text{ on } \Gamma ,

has exactly one solution; that solution is $C^2$ , and it is the function produced by steps 3 to 5 of the method. The regularity of the solution is the conclusion of Evans's Theorem 2, and the uniqueness clause is his specialisation of that theorem to the quasilinear case.

The proof below is Evans's route run for the three-equation system (2) alone, which is all the quasilinear case needs, and it uses less than his hypothesis: $C^1$ coefficients, a $C^1$ curve and $C^1$ data are enough for every step, and what they yield is a solution of class $C^1$ , unique among $C^1$ solutions on $V$ . The $C^2$ conclusion stated above is the one Evans proves, under his smoothness hypothesis.

For a linear equation the hypothesis simplifies, because $a$ and $b$ do not depend on $z$ : the condition is that the vector field $(a,b)$ is nowhere tangent to $\Gamma$ , and it involves neither the data $g$ nor the derivatives of $u$ . In that case Stakgold and Holst state, for the linear version of their first-order equation and without proof, a stronger conclusion: when $\Gamma$ is nowhere characteristic, a unique solution exists on the whole characteristic "strip" — their word, in quotation marks, for the region of the plane that the two characteristics through the ends of $\Gamma$ fence off between them — and not merely near one point, and it depends continuously on the data there. They add that a jump in the prescribed values of $u$ at a point of $\Gamma$ travels along the characteristic through that point, so that the solution is then a weak one in the sense of distributions.

Idea of the proofOriginal

Equation (1) says that one particular directional derivative of $u$ , in the direction $(a,b)$ , equals $c$ . So along a curve whose tangent is that direction, $u$ satisfies an ordinary differential equation. The three equations (2) are exactly that curve and that ordinary differential equation written together, and they close because the direction depends only on $x, y$ and the value of $u$ , never on its derivatives. Two things then have to be checked. First, that the characteristics issuing from $\Gamma$ fill out a neighbourhood without overlapping, which is the inverse function theorem applied to $(s,\tau) \mapsto (x,y)$ , and whose Jacobian determinant is precisely the non-characteristic quantity. Second, that the function so assembled really solves the equation: one line of the chain rule gives this, because $\partial z/\partial s$ is both $c$ and the directional derivative $a u_x + b u_y$ .

ProofAdaptedAdapted from Evans (2010), §3.2.4, Lemma 2 and Thm 2, pp. 102–105

Step 1 (the characteristic strips exist and are smooth in both variables). The right-hand sides of (2) are $C^1$ , hence locally Lipschitz in $(x,y,z)$ , so by the existence and uniqueness theorem for ordinary differential equations there is, for each $\tau$ near $\tau_0$ , a unique solution $s \mapsto \bigl(x(s,\tau), y(s,\tau), z(s,\tau)\bigr)$ of (2) with the initial values of step 3 of the method, defined for $|s| < s_0$ with $s_0 > 0$ independent of $\tau$ on a small enough $\tau$ -interval. Because the initial data $f, h, g$ are $C^1$ in $\tau$ and the right-hand sides are $C^1$ , the solution depends $C^1$ on $\tau$ as well.

Step 2 (the map $(s,\tau) \mapsto (x,y)$ is locally invertible). At $s = 0$ we have $x(0,\tau) = f(\tau)$ and $y(0,\tau) = h(\tau)$ , so $x_\tau = f'$ , $y_\tau = h'$ there, while (2) gives $x_s = a$ , $y_s = b$ evaluated on the data. Hence

\det \frac{\partial(x,y)}{\partial(s,\tau)}\bigg|_{(0,\tau_0)} = x_s y_\tau - x_\tau y_s = a\,h'(\tau_0) - b\,f'(\tau_0) ,

which is nonzero by the non-characteristic hypothesis. By the inverse function theorem there are a neighbourhood $W$ of $(0,\tau_0)$ and a neighbourhood $V$ of $P_0$ such that $(s,\tau) \mapsto (x,y)$ is a $C^1$ bijection of $W$ onto $V$ with $C^1$ inverse $(x,y) \mapsto (s(x,y), \tau(x,y))$ .

Step 3 (the assembled function solves the equation). Define $u(x,y) = z(s(x,y),\tau(x,y))$ on $V$ ; it is $C^1$ as a composition of $C^1$ maps. Fix $(s,\tau) \in W$ and differentiate the identity $z(s,\tau) = u\bigl(x(s,\tau), y(s,\tau)\bigr)$ with respect to $s$ . By the chain rule and then by (2),

\frac{\partial z}{\partial s} = u_x\,\frac{\partial x}{\partial s} + u_y\,\frac{\partial y}{\partial s} = a\,u_x + b\,u_y ,

all three functions $a, b$ and the derivatives being evaluated at the point $(x(s,\tau), y(s,\tau))$ with $z = u$ there. But the third equation of (2) says $\partial z/\partial s = c$ . Hence $a u_x + b u_y = c$ at every point of $V$ , which is (1). Putting $s = 0$ gives $u(f(\tau), h(\tau)) = z(0,\tau) = g(\tau)$ , so the data are matched.

Step 4 (uniqueness). Let $\tilde u$ be any $C^1$ solution on $V$ with $\tilde u = g$ on $\Gamma$ . Fix $\tau$ and solve the ordinary differential equations

\frac{d\tilde x}{ds} = a\bigl(\tilde x, \tilde y, \tilde u(\tilde x, \tilde y)\bigr), \qquad \frac{d\tilde y}{ds} = b\bigl(\tilde x, \tilde y, \tilde u(\tilde x, \tilde y)\bigr),

with $\tilde x(0) = f(\tau)$ , $\tilde y(0) = h(\tau)$ ; the right-hand sides are $C^1$ in $(\tilde x, \tilde y)$ because $a, b$ and $\tilde u$ are, so a unique solution exists for small $|s|$ . Put $\tilde z(s) = \tilde u(\tilde x(s), \tilde y(s))$ . By the chain rule and then because $\tilde u$ satisfies (1),

\frac{d\tilde z}{ds} = a\,\tilde u_x + b\,\tilde u_y = c\bigl(\tilde x, \tilde y, \tilde z\bigr) .

So $(\tilde x, \tilde y, \tilde z)$ solves the system (2) with the same initial values as $(x(\cdot,\tau), y(\cdot,\tau), z(\cdot,\tau))$ . The right-hand sides of (2) are locally Lipschitz, so the two solutions coincide. In particular $\tilde u(x(s,\tau), y(s,\tau)) = z(s,\tau) = u(x(s,\tau),y(s,\tau))$ , and since by Step 2 these points fill $V$ , $\tilde u = u$ on $V$ . $\square$

RemarkOriginal

Where the hypotheses are used, and what fails without them.

The non-characteristic condition is used only in Step 2, and nowhere else; it is exactly what makes the Jacobian nonzero. If it fails at $\tau_0$ , the characteristic through that point of $\Gamma$ is tangent to $\Gamma$ , the map $(s,\tau) \mapsto (x,y)$ is singular, and no neighbourhood of $P_0$ is swept out. Arfken draws the three consequences explicitly: data on a curve that is a characteristic are either inconsistent with the equation or merely redundant, and in neither case do they determine $u$ off that curve; a data curve meeting the same characteristic twice generally over-determines the problem, so no solution exists; and only a non-characteristic curve lets a boundary condition reach points off it.
Smoothness of the coefficients. Evans assumes it throughout, and he needs it for his $C^2$ conclusion. The proof above needs only that $a, b, c$ are $C^1$ , and it uses that twice, both times to get local Lipschitz continuity so that the ordinary differential equations have unique solutions — in Step 1 for existence and in Step 4 for uniqueness. Continuity alone would give existence but not uniqueness, and then the construction would not single out a solution.
The conclusion is local. This is not laziness. For a linear equation the characteristics are fixed curves, determined by $(a,b)$ alone, and the solution extends along each of them as far as they go, which is Stakgold and Holst's characteristic "strip"; the obstruction is then geometric, arising where the vector field $(a,b)$ has a zero inside the region or where a characteristic re-enters through the data curve. For a genuinely quasilinear equation the characteristics depend on the data through $z$ , and projected characteristics starting at different points of $\Gamma$ can intersect at a finite distance. Where two of them carrying different values of $u$ meet, no single-valued $C^1$ solution can exist: the derivatives of $u$ blow up as that time is approached, and the solution must be continued in a weaker sense (see conservation laws and shocks).
No smoothness of the data is needed for the formula. Steps 3 to 5 produce a function whenever $g$ is merely continuous, or even discontinuous, and it is the only sensible candidate; it is then a weak, not a classical, solution. Stakgold and Holst record that a jump in $g$ at a point of $\Gamma$ travels along the characteristic through that point and nowhere else.

Example(The transport equation)AdaptedAdapted from Evans (2010), §2.1, Eqs. (2)–(5), pp. 18–19

Let $b \in \mathbb{R}^n$ be a constant vector and solve

u_t + b\cdot Du = f \ \text{ in } \mathbb{R}^n\times(0,\infty), \qquad u = g \ \text{ on } \mathbb{R}^n\times\{t=0\} ,

where $Du$ is the gradient in $x$ alone. In the notation of (1), with $(x,t)$ in place of $(x,y)$ , the coefficient vector of the derivatives is $(b, 1) \in \mathbb{R}^{n+1}$ and the right-hand side is $f$ .

Data curve and the non-characteristic check. The data surface is $\Gamma = \{t = 0\}$ , with unit normal $(0,\dots,0,1)$ . The coefficient vector has last component $1$ , so $(b,1)\cdot(0,1) = 1 \ne 0$ everywhere on $\Gamma$ : the condition of step 2 holds for every $b$ and every point, which is why the initial-value problem for transport is always solvable.
Characteristics. The characteristic equations are $\dot x = b$ , $\dot t = 1$ , $\dot z = f$ . The first two give $t = s$ and $x = x^0 + sb$ , so we may use $t$ as the parameter: the projected characteristics are the straight lines through $(x^0, 0)$ with direction $(b,1)$ .
The equation along a characteristic. On the line through $(x,t)$ put $\zeta(s) = u(x + (s-t)b,\ s)$ , which traces the values of $u$ along that line and equals $u(x,t)$ at $s = t$ and $g(x - tb)$ at $s = 0$ . By the chain rule $\zeta'(s) = u_t + b\cdot Du$ evaluated on the line, which is $f(x + (s-t)b,\ s)$ .
Integrate. Integrating $\zeta'$ from $0$ to $t$ ,
$u(x,t) = g(x - tb) + \int_0^t f\bigl(x + (s-t)b,\ s\bigr)\,ds .$
For $f \equiv 0$ this is $u(x,t) = g(x - tb)$ : the initial profile is carried rigidly at velocity $b$ , unchanged in shape.
A number. Take $n = 2$ , $b = (3,-1)$ , $f \equiv 0$ and $g(x) = e^{-|x|^2}$ . Then $u(x,t) = e^{-|x - tb|^2}$ , and at $x = (1,2)$ , $t = 1$ we get $x - b = (-2, 3)$ , so $u = e^{-13} = 2.2603 \times 10^{-6}$ , while $u$ at the same point at $t = 0$ is $e^{-5} = 6.7379\times 10^{-3}$ : the pulse has moved away.
Regularity. The formula with $f \equiv 0$ defines a function for any $g$ whatever. It is a classical solution exactly when $g \in C^1$ ; for rougher $g$ it is the weak solution, and a jump in $g$ at $x^0$ travels along the line $x = x^0 + tb$ .

ExampleOriginal

(a) A linear equation with absorption and a source. Solve

u_t + 2u_x + 3u = e^{-x}, \qquad u(x,0) = \cos x .

Form and check. Here $(a,b) = (2,1)$ in the variables $(x,t)$ , and $c = e^{-x} - 3u$ . The data curve $t = 0$ has tangent $(1,0)$ , so $J = 2\cdot 0 - 1\cdot 1 = -1 \ne 0$ : non-characteristic everywhere.
Characteristics. $\dot x = 2$ , $\dot t = 1$ give $t = s$ , $x = x_0 + 2t$ , so the projected characteristics are the lines $x - 2t = x_0$ .
The ordinary differential equation along a characteristic. With $z(t) = u(x_0 + 2t, t)$ ,
$\frac{dz}{dt} = u_t + 2u_x = e^{-x} - 3u = e^{-x_0 - 2t} - 3z, \qquad z(0) = \cos x_0 .$
Solve it. Multiplying by the integrating factor $e^{3t}$ gives $(e^{3t}z)' = e^{t - x_0}$ , so $e^{3t}z = \cos x_0 + e^{-x_0}(e^{t}-1)$ and
$z(t) = e^{-3t}\cos x_0 + e^{-x_0}\bigl(e^{-2t} - e^{-3t}\bigr) .$
Return to $(x,t)$ . Substituting $x_0 = x - 2t$ and simplifying $e^{-(x-2t)}(e^{-2t}-e^{-3t}) = e^{-x}(1 - e^{-t})$ ,
$u(x,t) = e^{-3t}\cos(x - 2t) + e^{-x}\bigl(1 - e^{-t}\bigr) .$
Check. $u(x,0) = \cos x$ , and substituting into the equation gives $u_t + 2u_x + 3u = e^{-x}$ identically. As $t \to \infty$ the transient $e^{-3t}\cos(x-2t)$ dies and $u \to e^{-x}$ , the steady state of $2u_x + 3u = e^{-x}$ . At $x = 1$ , $t = \tfrac12$ the two terms are $e^{-3/2}\cos 0 = e^{-3/2}$ and $e^{-1}(1 - e^{-1/2}) = e^{-1} - e^{-3/2}$ , which cancel to give exactly $u(1,\tfrac12) = e^{-1} = 0.3678794$ .

(b) A genuinely quasilinear equation. Solve $u_t + u\,u_x = 0$ with $u(x,0) = 1 + x$ .

Characteristics. $\dot t = 1$ , $\dot x = z$ , $\dot z = 0$ . So $z \equiv 1 + x_0$ is constant along each characteristic, and then $x = x_0 + (1+x_0)t$ . The non-characteristic quantity on $t = 0$ is again $1 \ne 0$ .
Invert. From $x = x_0(1+t) + t$ we get $x_0 = (x-t)/(1+t)$ , which is legitimate for $t > -1$ .
Solution. $u = 1 + x_0 = \dfrac{1+x}{1+t}$ . Substituting back, $u_t = -(1+x)/(1+t)^2$ and $u\,u_x = (1+x)/(1+t)^2$ , so the equation holds, and $u(x,0) = 1+x$ .
Where it holds, and why. The characteristics are lines of slope $dx/dt = 1 + x_0$ , which increases with $x_0$ : they fan out for $t > 0$ and never meet, so the formula is valid for all $t > 0$ and the solution flattens like $1/(1+t)$ . Extended backwards they converge, and all of them pass through $(x,t) = (-1,-1)$ ; the solution indeed blows up as $t \to -1^+$ . Had the initial data been decreasing instead, the characteristics would have crossed at a positive time and the classical solution would have ceased to exist there.

(c) Data on a characteristic. For $2u_x + 3u_y = 0$ the characteristics are the lines $3x - 2y = \text{const}$ . Prescribing $u$ on the line $\Gamma = \{3x = 2y\}$ , parametrised by $(2\tau, 3\tau)$ , gives $J = a h' - b f' = 2\cdot 3 - 3\cdot 2 = 0$ : the data curve is characteristic at every point. The equation forces $u$ to be constant along $\Gamma$ , so prescribing $u = \tau$ there admits no solution at all, while prescribing the constant $u = 5$ admits infinitely many, namely $u = F(3x-2y)$ for every $C^1$ function $F$ with $F(0) = 5$ . The two failure modes — no existence and no uniqueness — appear in the same example, depending only on the data.

ApplicationOriginal

The transport equation is the advection step of every atmosphere and ocean model. A tracer — water vapour, ozone, a pollutant, a dye — of concentration $u$ carried by a prescribed wind field $\mathbf{v}(x,t)$ satisfies $u_t + \mathbf{v}\cdot Du = S$ , with $S$ the sources and sinks; its characteristics are the trajectories of air parcels, and solving along them is what "Lagrangian" and "semi-Lagrangian" advection schemes do literally: to find $u$ at a grid point at the new time level, they trace the characteristic backwards to the previous level and interpolate there. The method also explains the Courant condition for an explicit grid-based scheme, which says that the numerical domain of dependence of a grid point must contain the characteristic through it; violating it means computing a value from data that the exact solution does not depend on, and the result is unstable (see von Neumann stability analysis).

Elsewhere, the eikonal picture of geometrical optics and the ray equations of seismology are characteristic systems of first-order equations; the Hamilton-Jacobi equation of mechanics is one too, with the characteristics being the trajectories in phase space (see the Hamilton-Jacobi equation). In traffic and in gas dynamics the quasilinear case is the physical one, and part (b) of the example is the whole story in miniature: where characteristics spread, a smooth solution persists and flattens; where they converge, the gradient steepens without bound and a shock forms.

Method · advanced

Conservation laws and shocks

Adapted

Add to PDF

Shows how a scalar conservation law in one space dimension loses its classical solution at a computable breaking time, gives the weak formulation that lets the solution be continued past it, derives and proves the Rankine-Hugoniot speed of the resulting jump and the entropy condition that picks out the physical one, and carries a Burgers ramp through from the data to a checked shock position.

Prerequisites: First-order equations and characteristics, Implicit and inverse function theorems

Definition(Scalar conservation laws, flux, and Burgers' equation)AdaptedAdapted from Evans (2010), §3.4, Eq. (1), p. 130, and Example 1, Eq. (13), p. 134; Strauss (2008), §14.1, Eq. (1), p. 380, the three lettered principles $(\alpha)$, $(\beta)$, $(\gamma)$ of Example 3, p. 382, Eq. (13), p. 384, and the general-flux slope with Eq. (14), p. 385

A scalar conservation law in one space dimension is an equation in divergence form

u_t + \bigl(F(u)\bigr)_x = 0 , \qquad x \in \mathbb{R}, \ t > 0 , \tag{1}

for an unknown $u = u(x,t)$ and a given flux function $F : \mathbb{R} \to \mathbb{R}$ . We always pair (1) with initial data

u(x,0) = g(x) , \qquad x \in \mathbb{R} . \tag{2}

Equation (1) is Evans's Eq. (1); the name records the physics, since $u$ is a density whose total amount changes only through the flux $F(u)$ at the ends of an interval.

Wherever $u$ is $C^1$ the chain rule turns (1) into the quasilinear equation

u_t + a(u)\,u_x = 0 , \qquad a := F' , \tag{3}

which is Strauss's Eq. (1). We call $a(u)$ the characteristic speed. The two forms are equivalent for smooth $u$ and they are not equivalent for discontinuous $u$ : only (1) survives the passage to weak solutions below, because only (1) says what is conserved. This is why the flux, and not just the speed, has to be named.

Burgers' equation is the case $F(u) = u^2/2$ , so $a(u) = u$ :

u_t + \Bigl(\tfrac{1}{2}u^2\Bigr)_x = 0 \quad\text{equivalently}\quad u_t + u\,u_x = 0 . \tag{4}

Following Evans we call $F$ uniformly convex when $F'' \ge \theta$ for some constant $\theta > 0$ ; Burgers' flux has $F'' = 1$ , so it qualifies. Under uniform convexity $a = F'$ is strictly increasing, a fact used repeatedly below.

Equation (3) is the quasilinear equation of first-order equations and characteristics with $(a,b,c) = (a(u), 1, 0)$ in the variables $(x,t)$ . Its characteristic equations are therefore $\dot x = a(z)$ , $\dot t = 1$ , $\dot z = 0$ , and from them come the three facts that organise everything here: $u$ is constant along each projected characteristic, each projected characteristic is consequently a straight line, and the slope $dx/dt$ of that line is the characteristic speed $a(u)$ carried on it. Strauss sets out the three facts as the lettered principles $(\alpha)$ , $(\beta)$ , $(\gamma)$ on p. 382, but he derives them there inside his Example 3, which is Burgers' equation, $a(u) = u$ ; his third principle therefore equates the slope with $u$ itself, an identification special to $a = \mathrm{id}$ and not available for a general flux. The general-flux version is his p. 385, where the line through $(x,t)$ and $(z,0)$ has slope $(x-z)/t = a(u(x,t)) = a(\varphi(z))$ , and it is the version used throughout this unit. The line through $(y,0)$ is

x = y + t\,a\bigl(g(y)\bigr) , \qquad u = g(y) \quad\text{on it.} \tag{5}

Inverting (5) for $y$ gives $u(x,t) = g(y(x,t))$ , which is Strauss's implicit solution (14). The data curve $\{t=0\}$ is never characteristic, since the coefficient of $u_t$ is $1 \ne 0$ , so the only obstruction to this construction is that two lines (5) may meet.

Method(Breaking: how long the classical solution lasts)AdaptedAdapted from Strauss (2008), §14.1, Eq. (14) and the crossing criterion $a(\varphi(z)) \le a(\varphi(w))$ for $z \le w$, pp. 383–385

Assume $F \in C^2(\mathbb{R})$ and $g \in C^1(\mathbb{R})$ with $g$ , $g'$ , $a \circ g$ and $(a\circ g)'$ all bounded, where $a = F'$ . Put

m := \inf_{y \in \mathbb{R}} \, (a \circ g)'(y) = \inf_{y\in\mathbb{R}} a'\bigl(g(y)\bigr) g'(y) .

No crossing. If $m \ge 0$ then for each $t > 0$ the map $y \mapsto y + t\,a(g(y))$ is a strictly increasing bijection of $\mathbb{R}$ onto itself, and (5) defines a bounded $C^1$ solution of (1)–(2) on the whole half-plane $t > 0$ . Strauss's criterion for this, that $a(g(\cdot))$ be nondecreasing, is exactly $m \ge 0$ ; he calls the resulting solution an expansive or rarefaction wave.
Breaking. If $m < 0$ , set the breaking time
$t_* := -\frac{1}{m} \;=\; \frac{1}{\,\bigl|\inf_y (a\circ g)'(y)\bigr|\,} \;>\; 0 .$
Then (5) defines a bounded $C^1$ solution on the strip $0 \le t < t_*$ , and it is the only bounded $C^1$ solution there. If the infimum is attained at some $y_*$ , then
$\sup_{x \in \mathbb{R}} |u_x(x,t)| \longrightarrow \infty \qquad\text{as } t \uparrow t_* ,$
so no bounded $C^1$ solution exists on any strip $0 \le t < T$ with $T > t_*$ , while $\sup_x|u(x,t)| = \sup_x |g|$ for every $t < t_*$ .

The quantity that blows up is the slope, not the solution. The derivation of the explicit number $t_*$ from Strauss's implicit formula is carried out in the proof below; Strauss states only the sign criterion that separates case 1 from case 2 and describes the breaking in case 2 geometrically, as the crest of the wave overtaking its foot.

Idea of the proofOriginal

Read (5) backwards. At time $t$ the characteristic that started at $y$ has arrived at $x = y + t\,a(g(y))$ , so the whole initial line has been stretched by the map $y \mapsto x$ , and the stretching factor is $1 + t\,(a\circ g)'(y)$ . The solution is obtained by undoing that map, which is possible exactly while the factor stays positive. Since $u$ keeps the value $g(y)$ along the characteristic, its slope is the initial slope divided by the stretching factor; where the factor shrinks to zero the graph of $u$ has been compressed to a vertical tangent. The first time this happens is the first zero of $1 + t\,(a\circ g)'(y)$ over all $y$ , and that is $-1/m$ .

ProofAdaptedAdapted from Strauss (2008), §14.1, pp. 383–385

Write $\alpha := a \circ g$ , which is $C^1$ and bounded with bounded derivative, and define

\Phi(y,t) := y + t\,\alpha(y) .

Step 1 (the stretching factor). $\partial_y \Phi(y,t) = 1 + t\,\alpha'(y)$ . If $m \ge 0$ this is at least $1$ for all $t > 0$ ; if $m < 0$ and $0 \le t < t_* = -1/m$ it is at least $1 + t\,m > 0$ . In either case, on the stated range of $t$ the function $y \mapsto \Phi(y,t)$ is strictly increasing.

Step 2 (it is onto). $\alpha$ is bounded, say $|\alpha| \le M$ , so $\Phi(y,t) \ge y - tM \to \infty$ and $\Phi(y,t) \le y + tM \to -\infty$ as $y \to \pm\infty$ . A continuous strictly increasing function with these limits is a bijection of $\mathbb{R}$ onto $\mathbb{R}$ . Hence for each $t$ in range and each $x$ there is exactly one root $y = y(x,t)$ of $\Phi(y,t) = x$ .

Step 3 (the root is a $C^1$ function of $(x,t)$ ). The map $(y,t) \mapsto (\Phi(y,t), t)$ is $C^1$ with Jacobian determinant $\partial_y\Phi = 1 + t\alpha'(y) \ne 0$ at every point of the strip. By the inverse function theorem it has a $C^1$ inverse near every point, and by Step 2 the inverse is globally single-valued; so $y(x,t)$ is $C^1$ on the strip.

Step 4 (the formula solves the problem). Put $u(x,t) := g(y(x,t))$ , a $C^1$ function bounded by $\sup|g|$ . Differentiating the identity $g(y) = u(\Phi(y,t), t)$ in $y$ and then in $t$ , at fixed $t$ and fixed $y$ respectively,

g'(y) = u_x \cdot \bigl(1 + t\alpha'(y)\bigr), \qquad 0 = u_x\,\alpha(y) + u_t .

The second identity is $u_t + a(g(y))\,u_x = 0$ , and $g(y) = u$ at that point, so $u_t + a(u)u_x = 0$ , which is (3) and hence (1) for a $C^1$ function. Setting $t = 0$ in $\Phi$ gives $y(x,0) = x$ and so $u(x,0) = g(x)$ . Every value taken by $u$ is a value of $g$ , which gives the bound $\sup_x|u(\cdot,t)| = \sup|g|$ .

Step 5 (the slope blows up at $t_*$ ). From the first identity of Step 4,

u_x\bigl(\Phi(y,t),t\bigr) = \frac{g'(y)}{1 + t\,\alpha'(y)} . \tag{6}

Suppose $m < 0$ is attained at $y_*$ , so $\alpha'(y_*) = a'(g(y_*))g'(y_*) = m$ . A product is negative only if both factors are nonzero, so $g'(y_*) \ne 0$ . Evaluating (6) at $y_*$ gives $|u_x| = |g'(y_*)|/(1 + t m)$ , whose denominator decreases to $0$ as $t \uparrow t_* = -1/m$ . Hence $\sup_x|u_x(\cdot,t)| \to \infty$ . A $C^1$ solution on a strip reaching beyond $t_*$ would, by the uniqueness of Step 6, agree with $u$ below $t_*$ and so have unbounded derivative on a compact set, which is impossible.

Step 6 (uniqueness). Let $\tilde u$ be a bounded $C^1$ solution on the strip. Fix a point $(x,t)$ in the strip and solve $\dot\xi(\tau) = a(\tilde u(\xi(\tau),\tau))$ backwards from $\xi(t) = x$ . The right-hand side is $C^1$ in $\xi$ , hence locally Lipschitz, so the solution is unique, and it is bounded because $\tilde u$ is bounded and $a$ is continuous; therefore it continues down to $\tau = 0$ . Along it, $\frac{d}{d\tau}\tilde u(\xi,\tau) = \tilde u_t + a(\tilde u)\tilde u_x = 0$ , so $\tilde u$ keeps the value $g(\xi(0))$ and $\dot\xi$ keeps the value $a(g(\xi(0)))$ . Integrating, $x = \xi(0) + t\,a(g(\xi(0)))$ , so $\xi(0)$ is a root of $\Phi(\cdot,t) = x$ ; by Step 2 it is the root $y(x,t)$ , and $\tilde u(x,t) = g(y(x,t)) = u(x,t)$ . This is the backward-characteristic argument of first-order equations and characteristics, run globally on the strip rather than near one point. $\square$

RemarkOriginal

Where the hypotheses are used, and what fails without them.

$F \in C^2$ and $g \in C^1$ are used in Step 3, to make $\alpha = F' \circ g$ continuously differentiable so that the inverse function theorem applies, and in Step 6, to make $a(\tilde u)$ locally Lipschitz so that the backward characteristic is unique. With $F$ merely $C^1$ the map $\Phi$ need not be differentiable in $y$ and the breaking time is not defined by the formula above.
Boundedness of $\alpha$ is used only in Step 2, and only to get surjectivity. Without it the characteristics may fail to cover the half-plane, so the formula would define a solution on a smaller region even before any crossing.
The infimum being attained is used only in Step 5. If it is approached but not attained the solution still ceases to be $C^1$ at $t_*$ , but the blow-up is then at spatial infinity rather than at a point.
$m < 0$ is not an accident. If $F$ is uniformly convex, $a' = F'' > 0$ , so $m < 0$ as soon as $g$ is strictly decreasing anywhere. Thus for Burgers' equation any initial profile with a descending part breaks in finite time, no matter how smooth or how small. This is the characteristic difference from the linear transport equation, where $a$ is constant, $\alpha' \equiv 0$ and $t_* = \infty$ .
What the theorem does not say. It says nothing about $t \ge t_*$ . The solution does not become multivalued — a function cannot — but the construction (5) does: past $t_*$ it assigns several values to the same point, which is the folded profile Strauss draws as a breaking wave. Continuing past $t_*$ therefore needs a notion of solution that does not mention $u_x$ at all.

Definition(Integral solutions)AdaptedAdapted from Evans (2010), §3.4.1a, Eqs. (2)–(4) and the definition of an integral solution, pp. 131–132; Strauss (2008), §14.1, Eq. (15), p. 386

Call $v : \mathbb{R}\times[0,\infty) \to \mathbb{R}$ a test function if it is smooth and vanishes outside a bounded set. Multiply (1) by $v$ , integrate over the half-plane and move both derivatives onto $v$ : the $x$ -integration by parts has no boundary term because $v$ has compact support, and the $t$ -integration by parts leaves the boundary term at $t = 0$ , where $u = g$ . The result is an identity that contains no derivative of $u$ , and this is taken as the definition.

A function $u \in L^\infty(\mathbb{R}\times(0,\infty))$ is an integral solution (a weak solution) of (1)–(2) if

\int_0^\infty\!\!\int_{-\infty}^{\infty} \bigl( u\,v_t + F(u)\,v_x \bigr)\,dx\,dt \;+\; \int_{-\infty}^{\infty} g(x)\,v(x,0)\,dx \;=\; 0 \tag{7}

for every test function $v$ . This is Evans's Eq. (4) and his definition. Note where the flux enters: it is $F(u)$ , not $a(u)u_x$ , that is differentiated, so (7) makes sense for a merely bounded $u$ , and a smooth integral solution is a classical solution — reverse the integrations by parts and use that $v$ is arbitrary.

Two conventions. Strauss's Eq. (15) is (7) without the data term: his test functions live in the open half-plane, so his identity constrains $u$ only for $t>0$ and the initial condition is imposed separately. Evans's (7) carries the initial data inside the identity. We follow Evans, so that "integral solution" below always means a solution of the initial-value problem; when only the first integral in (7) vanishes we say $u$ solves (1) weakly in the open half-plane. The two agree for a $u$ that attains $g$ continuously as $t \downarrow 0$ , and the results quoted from Strauss are unaffected, because the jump conditions he derives are local in $t > 0$ .

Method(The Rankine-Hugoniot condition)AdaptedAdapted from Evans (2010), §3.4.1a, Eqs. (10)–(12) with the jump notation, pp. 132–133; Strauss (2008), §14.1, Eq. (17), p. 387

Let $V \subset \mathbb{R}\times(0,\infty)$ be open and let $C = \{(x,t) : x = s(t)\}$ be a $C^1$ curve dividing $V$ into a left part $V_l$ and a right part $V_r$ . Suppose $u$ is a bounded integral solution of (1) which, together with its first derivatives, is uniformly continuous on $V_l$ and on $V_r$ separately. Write

u_l(t) := \lim_{x \uparrow s(t)} u(x,t), \qquad u_r(t) := \lim_{x \downarrow s(t)} u(x,t),

and use the jump brackets $[\,u\,] := u_l - u_r$ , $[\,F(u)\,] := F(u_l) - F(u_r)$ , and $\sigma := \dot s$ for the speed of the curve. Then

$u$ is a classical solution of (1) in $V_l$ and in $V_r$ ; and
at every point of $C$ ,
$[\,F(u)\,] \;=\; \sigma\,[\,u\,] , \qquad\text{equivalently}\qquad \sigma \;=\; \frac{F(u_l) - F(u_r)}{u_l - u_r} \ \ \text{ when } u_l \ne u_r . \tag{8}$

This is the Rankine-Hugoniot condition. It is Evans's Eq. (12); Strauss's Eq. (17) is the same equation written as a difference quotient, with his $u^-$ and $u^+$ for the limits from the left and from the right, so that his $u^- = u_l$ and $u^+ = u_r$ and his $s(t)$ is our $\sigma$ . The three quantities $u_l$ , $u_r$ and $\sigma$ may each vary along $C$ ; what (8) fixes is the relation between them at each point. Nothing in (8) requires $F$ to be convex.

Note what the condition is a statement about: the flux, not the characteristic speed. Replacing $F$ by a different flux with the same $F'$ on the two states would change nothing in (3) but would change $\sigma$ — one more reason the conservation form (1) is the primitive object.

Idea of the proofOriginal

A weak solution satisfies the conservation law in its book-keeping form: for any fixed interval, the amount of $u$ inside changes at a rate equal to the flux in at the left end minus the flux out at the right end. Freeze an interval straddling the jump. The rate of change splits into two pieces: the smooth evolution on each side of the jump, which is accounted for by the fluxes at the two ends and at the jump, and the sweeping of the discontinuity across the interval, which transfers an amount $\sigma\,[\,u\,]$ per unit time from the right state to the left state. For the book-keeping to balance, that transfer must equal the mismatch $[\,F(u)\,]$ in the flux across the jump. Writing that out is the whole proof.

ProofAdaptedAdapted from Evans (2010), §3.4.1a, pp. 131–133

Step 1 (the equation holds classically off $C$ ). Let $v$ be a test function supported in $V_l$ . Then only the first integral in (7) survives, and since $u$ is $C^1$ there we may integrate by parts back again with no boundary terms, getting $\int\!\!\int_{V_l} (u_t + F(u)_x) v \,dx\,dt = 0$ . The integrand's first factor is continuous on $V_l$ and $v$ is arbitrary, so $u_t + F(u)_x = 0$ throughout $V_l$ . The same argument applies in $V_r$ . This is part 1.

Step 2 (the control-volume identity). Fix $a < b$ and $0 < t_1 < t_2$ with the closed rectangle $R = [a,b]\times[t_1,t_2]$ contained in $V$ , and with the vertical sides $x = a$ and $x = b$ not meeting $C$ . We claim

\int_a^b u(x,t_2)\,dx - \int_a^b u(x,t_1)\,dx = \int_{t_1}^{t_2} \Bigl[\, F\bigl(u(a,t)\bigr) - F\bigl(u(b,t)\bigr) \,\Bigr]\,dt . \tag{9}

To see it, take $v(x,t) = \chi_\varepsilon(x)\,\theta_\varepsilon(t)$ in (7), where $\chi_\varepsilon$ is smooth, equal to $1$ on $[a+\varepsilon, b-\varepsilon]$ , equal to $0$ outside $[a,b]$ , and monotone on each of the two transition intervals, and $\theta_\varepsilon$ is built the same way from $[t_1,t_2]$ . Since $R \subset V$ and $t_1>0$ , this $v$ is a test function and the data term in (7) is absent. Now $\theta_\varepsilon'$ integrates to $+1$ over the first transition and to $-1$ over the second, and $u$ is bounded and continuous in $t$ at each fixed $x$ off $C$ ; so the first double integral in (7) tends to $\int_a^b u(x,t_1)dx - \int_a^b u(x,t_2)dx$ . Likewise $\chi_\varepsilon'$ concentrates at $a$ with total mass $+1$ and at $b$ with total mass $-1$ , and $x \mapsto F(u(x,t))$ is continuous near $x = a$ and $x = b$ because those lines avoid $C$ ; so the second double integral tends to $\int_{t_1}^{t_2}[F(u(a,t)) - F(u(b,t))]\,dt$ . Both passages are dominated by the bound on $u$ and the continuity of $F$ on a bounded set. Adding gives (9).

Step 3 (differentiate in time). Fix $t$ with $t_1 < t < t_2$ . Dividing (9) by $t_2 - t_1$ and letting $t_2 \downarrow t_1 = t$ , both sides have limits and

\frac{d}{dt}\int_a^b u(x,t)\,dx \;=\; F\bigl(u(a,t)\bigr) - F\bigl(u(b,t)\bigr) . \tag{10}

Step 4 (split the interval at the jump). Assume $a < s(t) < b$ . Then

\int_a^b u\,dx = \int_a^{s(t)} u\,dx + \int_{s(t)}^{b} u\,dx ,

and each piece may be differentiated by the Leibniz rule, legitimately because $u$ and $u_t$ extend continuously to the closure of each side and $s$ is $C^1$ :

\frac{d}{dt}\int_a^b u\,dx = \sigma\,u_l - \sigma\,u_r + \int_a^{s(t)} u_t\,dx + \int_{s(t)}^{b} u_t\,dx .

By Step 1, $u_t = -F(u)_x$ on each side separately, so the two integrals evaluate to $-[F(u_l) - F(u(a,t))]$ and $-[F(u(b,t)) - F(u_r)]$ . Hence

\frac{d}{dt}\int_a^b u\,dx = \sigma\bigl(u_l - u_r\bigr) - F(u_l) + F(u_r) + F\bigl(u(a,t)\bigr) - F\bigl(u(b,t)\bigr) .

Step 5 (compare). Subtracting (10) from this leaves $0 = \sigma(u_l - u_r) - F(u_l) + F(u_r)$ , which is (8). The point $t$ was arbitrary. $\square$

Method(The entropy condition and the definition of a shock)AdaptedAdapted from Evans (2010), §3.4.1b, Eqs. (16)–(19), pp. 135–137, and §3.4.3b, Definition and Thm 3, p. 143; Strauss (2008), §14.1, Eq. (18), pp. 388–389

Integral solutions are not unique: for one two-state problem for Burgers' equation Evans exhibits two, a discontinuous one whose jump he checks against (8) and a continuous one with no jump at all. A further condition is therefore needed to single out the physical solution.

Entropy inequalities. Along a curve of discontinuity of a piecewise-smooth integral solution, Evans requires

F'(u_l) \;>\; \sigma \;>\; F'(u_r) . \tag{11}

These are the entropy inequalities, so named by analogy with the thermodynamic requirement that entropy not decrease. A curve of discontinuity is called a shock when both (8) and (11) hold on it. Strauss's criterion (18) is the same pair of inequalities written with $a = F'$ and his left and right limits, and he too declares the definition of a shock wave complete only once both conditions are imposed.

Why these inequalities. The characteristic through a point of $C$ carries the speed $F'(u_l)$ on the left and $F'(u_r)$ on the right. Condition (11) says both characteristics run into the curve as $t$ increases, so that going backwards in time from any point one never crosses a discontinuity; Evans obtains (11) exactly by imposing that backward-tracing property. Physically, as Strauss puts it, the wave behind the shock is overtaking the wave in front of it.

The convex case. If in addition $F$ is uniformly convex, then $F'$ is strictly increasing and (11) is equivalent to the single inequality

u_l \;>\; u_r \tag{12}

along the shock curve. For Burgers' equation ( $F' = \mathrm{id}$ ) condition (11) reads $u_l > \sigma > u_r$ , and (8) gives $\sigma = (u_l+u_r)/2$ , the average of the two states; the inequalities are then automatic once $u_l > u_r$ .

Uniqueness. Evans turns the pointwise inequality into a condition that a merely bounded function can satisfy: $u \in L^\infty$ is an entropy solution of (1)–(2) if it satisfies (7) for every test function and, for some constant $C \ge 0$ ,

u(x+z,t) - u(x,t) \;\le\; C\Bigl(1 + \frac{1}{t}\Bigr) z \qquad \text{for a.e. } x,z \in \mathbb{R},\ t>0,\ z>0 . \tag{13}

His Theorem 3 states that if $F$ is convex and smooth, then there is at most one entropy solution of (1)–(2), up to a set of measure zero. Condition (13) is a one-sided bound on how fast $u$ may increase in $x$ ; at a jump it forces $u_l \ge u_r$ , which is (12).

Riemann's problem. For $F$ uniformly convex and $C^2$ , with $G := (F')^{-1}$ , Evans's Theorem 4 gives the entropy solution for the two-state data $g = u_l$ for $x<0$ , $g = u_r$ for $x>0$ :

u_l > u_r : \qquad u(x,t) = \begin{cases} u_l, & x/t < \sigma, \\ u_r, & x/t > \sigma,\end{cases} \tag{14a}

u_l < u_r : \qquad u(x,t) = \begin{cases} u_l, & x/t < F'(u_l), \\ G(x/t), & F'(u_l) < x/t < F'(u_r), \\ u_r, & x/t > F'(u_r), \end{cases} \tag{14b}

with $\sigma$ given by (8). The first is a shock of constant speed, the second a rarefaction wave. The two formulas (14a) and (14b) are the building blocks out of which the solution past a breaking time is assembled.

RemarkOriginal

Where the hypotheses are used, and what fails without them.

Piecewise smoothness with uniformly continuous one-sided limits is used in Steps 1 and 4 of the Rankine-Hugoniot proof: Step 1 needs $u$ to be $C^1$ in the interior of each side, and Step 4 needs $u$ and $u_t$ to extend continuously up to $C$ from each side so that the Leibniz rule may be applied. Without it a bounded integral solution can have a much wilder discontinuity set and (8) has nothing to be about.
The vertical sides of the rectangle avoiding $C$ is used only in Step 2, to make $x\mapsto F(u(x,t))$ continuous where $\chi'_\varepsilon$ concentrates. Any $a,b$ with $a < s(t) < b$ on the time interval considered will do.
The divergence form is used throughout, and which divergence form matters. Multiplying $u_t + uu_x = 0$ by $2u$ turns it into $(u^2)_t + \bigl(\tfrac23 u^3\bigr)_x = 0$ , a conservation law for $w = u^2$ with flux $\tfrac23 w^{3/2}$ . The two laws have exactly the same $C^1$ solutions, but not the same jump speeds. For the states $u_l = 3$ , $u_r = 1$ the first gives $\sigma = \tfrac12(u_l+u_r) = 2$ , while the second gives $\bigl(\tfrac23\bigr)(27-1)/(9-1) = \tfrac{13}{6}$ . Only one of the two can be the physically correct front, and which one is settled by which quantity is actually conserved, not by the differential equation. A conservation law is specified by its flux.
Convexity is needed for (12) but not for (8) or (11). The general inequalities (11) are the ones to use for a nonconvex or concave flux; the traffic-flow model in the application below has $F'' < 0$ , and there (11) forces $u_l < u_r$ , the opposite of (12).
Why a non-entropic jump has to be rejected. This is Evans's Example 2 (p. 135) with different states: take Burgers' equation with $u_l = 1$ and $u_r = 3$ in place of his $0$ and $1$ ; the arithmetic below is ours. The jump $u = 1$ for $x < 2t$ , $u = 3$ for $x > 2t$ satisfies (8), since $(F(1)-F(3))/(1-3) = (\tfrac12 - \tfrac92)/(-2) = 2$ , and it is an integral solution. But $F'(u_l) = 1 < 2$ , so (11) fails at once. The entropy solution is instead the rarefaction of (14b): $u = 1$ for $x < t$ , $u = x/t$ for $t < x < 3t$ , $u = 3$ for $x > 3t$ . Both have the same initial data, so without (11) the initial-value problem would not determine the solution. Observe that the characteristics leave the rejected jump on both sides, so the state behind it is not determined by the data at all — that is precisely what the backward-tracing requirement forbids.

Method(Solving a scalar conservation law through a shock)AdaptedAdapted from Evans (2010), §3.4.1, pp. 131–138, and §3.4.4, Thm 4, pp. 147–148; Strauss (2008), §14.1, pp. 383–389

Given $u_t + F(u)_x = 0$ with $u(x,0) = g(x)$ :

Put the equation in divergence form and identify $F$ and $a = F'$ . If the equation is handed to you as $u_t + a(u)u_x = 0$ , choose $F$ with $F' = a$ ; the constant of integration does not matter, but the choice of which quantity is conserved does.
Draw the characteristics. Each is the line $x = y + t\,a(g(y))$ carrying the constant value $g(y)$ .
Find the breaking time. Compute $m = \inf_y (a\circ g)'(y)$ . If $m \ge 0$ the classical solution is global and step 4 onwards is unnecessary. Otherwise $t_* = -1/m$ , and the classical solution is given implicitly by $u = g(y)$ , $x = y + t\,a(g(y))$ for $t < t_*$ .
Insert the right wave at each jump of the data using (14a)–(14b): where the data jumps so that the characteristics converge, a shock; where they diverge, a rarefaction fan $u = G(x/t)$ centred at the jump, with $G = (F')^{-1}$ .
Start the shock. At $t = t_*$ the characteristics first cross; take the crossing point as the initial point of the shock curve $x = s(t)$ .
Propagate the shock. Read off $u_l(t)$ and $u_r(t)$ as the limits of the solution already constructed on either side of $s(t)$ , and integrate the ordinary differential equation
$\dot s = \frac{F(u_l) - F(u_r)}{u_l - u_r} , \qquad s(t_*) = \text{the crossing point} .$
When both states are constant this is a straight line; when one side is a fan it is a genuine ordinary differential equation in $s$ and $t$ .
Check the entropy inequalities $F'(u_l) > \dot s > F'(u_r)$ at each time. If they fail, the jump you inserted is not a shock and step 4 was done wrongly.
Stop and restart at interactions. The construction is valid until the shock meets another shock, the edge of a fan, or a new breaking point; at such a time recompute $u_l$ , $u_r$ and go back to step 6.
Verify. Check the control-volume identity (9) on an interval containing the shock: the amount of $u$ in $[a,b]$ must change at the rate $F(u(a,t)) - F(u(b,t))$ . This catches an error in $F$ , in $\sigma$ or in the solution on either side.

Example(A Burgers ramp: breaking, shock and a conserved count)Original

Solve Burgers' equation $u_t + \bigl(\tfrac12 u^2\bigr)_x = 0$ with

g(x) = \begin{cases} 3, & x \le 0,\\ 3 - 2x, & 0 \le x \le 1,\\ 1, & x \ge 1.\end{cases}

Here $F(u) = u^2/2$ and $a(u) = F'(u) = u$ , so $a \circ g = g$ . The data are continuous and piecewise $C^1$ but not $C^1$ , because of the corners at $x = 0$ and $x = 1$ , so the hypotheses of the breaking statement are not literally met. The corners cause no trouble: the characteristic speeds match across them, so the lines (5) still fill the strip without gaps, and the construction goes through on each piece separately. The price is that what comes out is $C^1$ only off the two characteristics $x = 3t$ and $x = 1+t$ carried by the corners, and merely Lipschitz across those two lines; we call it the piecewise-classical solution. The breaking time is unaffected, since the stretching factor $1 + t\,g'(y)$ is computed piece by piece, and no number below depends on the distinction.

1. Characteristics. The line from $(y,0)$ is $x = y + t\,g(y)$ , carrying the value $g(y)$ . So

for $y \le 0$ : $x = y + 3t$ , speed $3$ ;
for $0 \le y \le 1$ : $x = y + t(3-2y) = y(1-2t) + 3t$ , speed $3-2y$ ;
for $y \ge 1$ : $x = y + t$ , speed $1$ .

2. Breaking time. $(a\circ g)' = g'$ equals $0$ outside $(0,1)$ and $-2$ inside it, so $m = -2$ and

t_* = -\frac1m = \frac12 .

3. Where the crossing happens. At $t = \tfrac12$ the middle family collapses: $x = y(1-2t) + 3t$ becomes $x = \tfrac32$ for every $y \in [0,1]$ . All the characteristics carrying the values between $1$ and $3$ arrive at the single point $(x,t) = \bigl(\tfrac32, \tfrac12\bigr)$ .

4. The piecewise-classical solution for $t < \tfrac12$ . Inverting $x = y(1-2t)+3t$ gives $y = (x-3t)/(1-2t)$ and hence $u = 3 - 2y = (3-2x)/(1-2t)$ . Altogether

u(x,t) = \begin{cases} 3, & x \le 3t,\\[2pt] \dfrac{3-2x}{1-2t}, & 3t \le x \le 1+t,\\[4pt] 1, & x \ge 1+t,\end{cases} \qquad 0 \le t < \tfrac12 .

The two branch points agree with the ends of the middle family, and the formula is continuous across them: at $x = 3t$ it gives $(3-6t)/(1-2t) = 3$ and at $x = 1+t$ it gives $(1-2t)/(1-2t) = 1$ . A check at an interior point: at $x = 1$ , $t = \tfrac14$ the formula gives

u\bigl(1,\tfrac14\bigr) = \frac{3-2}{1-\tfrac12} = 2 ,

and indeed the characteristic through $(1,\tfrac14)$ started at $y = (1 - \tfrac34)/\tfrac12 = \tfrac12$ , where $g = 3 - 1 = 2$ . The slope at that point is $u_x = -2/(1-2t) = -4$ , which does run to $-\infty$ as $t \uparrow \tfrac12$ , as Step 5 of the proof predicts.

5. The shock. For $t > \tfrac12$ the states on the two sides of the crossing point are the constants $u_l = 3$ and $u_r = 1$ . By (8),

\sigma = \frac{F(3) - F(1)}{3-1} = \frac{\tfrac92 - \tfrac12}{2} = 2 ,

constant, so the shock curve is the straight line through $(\tfrac32,\tfrac12)$ of slope $dx/dt = 2$ :

s(t) = \tfrac32 + 2\bigl(t - \tfrac12\bigr) = 2t + \tfrac12 .

6. Entropy. $F'(u_l) = 3$ , $F'(u_r) = 1$ , and $3 > 2 > 1$ : the inequalities (11) hold, so this is a genuine shock. (Equivalently, $u_l = 3 > 1 = u_r$ , since $F$ is uniformly convex with $F''=1$ .) The solution for $t \ge \tfrac12$ is therefore

u(x,t) = \begin{cases} 3, & x < 2t + \tfrac12,\\ 1, & x > 2t + \tfrac12. \end{cases}

Nothing further happens: the characteristics on the left have speed $3 > 2$ and on the right speed $1 < 2$ , so both run into the shock for all later times and no new interaction occurs.

7. A checkable number. At $t = 3$ the shock sits at

s(3) = 2\cdot 3 + \tfrac12 = 6.5 .

Verify this against the conservation identity (9) on the fixed interval $[-10, 10]$ , which contains the shock for all $0 \le t \le 3$ (at $t=3$ it is at $6.5 < 10$ ) and whose ends sit in the constant states. Write $N(t) := \int_{-10}^{10} \bigl(u(x,t) - 1\bigr)\,dx$ ; since the constant $1$ contributes nothing to the right-hand side of (9) beyond a cancelling pair, (9) gives

N(t) - N(0) = \int_0^t \bigl[F(u(-10,\tau)) - F(u(10,\tau))\bigr]\,d\tau = t\bigl(F(3) - F(1)\bigr) = 4t .

Now compute both sides directly. At $t=0$ ,

N(0) = \int_{-10}^{0} 2\,dx + \int_0^1 \bigl(2 - 2x\bigr)\,dx + 0 = 20 + 1 = 21 .

At $t = 3$ the profile takes the value $u - 1 = 2$ on $[-10, 6.5)$ and $0$ beyond, so

N(3) = 2\cdot\bigl(6.5 - (-10)\bigr) = 33 .

The difference is $33 - 21 = 12 = 4 \cdot 3$ , exactly as required. Had the shock been placed anywhere other than $x = 6.5$ at $t=3$ , the count would have been wrong: the Rankine-Hugoniot speed is precisely the speed that keeps the books balanced.

ApplicationOriginal

Traffic. The Lighthill-Whitham-Richards model treats the density $\rho(x,t)$ of vehicles on a single lane as a conserved quantity whose flux is $F(\rho) = v_m\rho\,(1 - \rho/\rho_m)$ : drivers choose the speed $v_m(1-\rho/\rho_m)$ , full speed on an empty road and zero at the jam density $\rho_m$ . Here $F'' = -2v_m/\rho_m < 0$ , so the flux is concave, not convex, and the entropy inequalities must be used in the general form (11).

Take $v_m = 100\ \mathrm{km/h}$ and $\rho_m = 180$ vehicles per km, and let traffic of density $\rho_0 = 36$ per km run into a standing queue at $\rho = \rho_m$ behind a red light. The queue's tail is a jump with $\rho_l = 36$ upstream and $\rho_r = 180$ downstream, and (8) gives

\sigma = \frac{F(36) - F(180)}{36 - 180} = \frac{2880 - 0}{-144} = -20\ \mathrm{km/h} :

the back of the queue travels upstream at $20\ \mathrm{km/h}$ while every car in it moves forwards or stands still. The entropy check is $F'(36) = 60 > -20 > -100 = F'(180)$ , so this is a shock. Over three minutes the tail retreats $1\ \mathrm{km}$ , so the queue grows by $180$ vehicles; and indeed $144$ cars arrive in that time at $2880$ per hour, while $36$ were already on that kilometre of road. The two counts agree, which is the control-volume check of step 9 of the method. When the light turns green the jump runs the other way, $\rho_l = \rho_m > \rho_r$ , the characteristics diverge, and (14b) supplies a rarefaction fan instead: this is why a queue dissolves gradually from the front while it builds abruptly at the back.

Gas dynamics. The equations of compressible flow are a system of conservation laws for mass, momentum and energy, and the jump relations across a shock front — the original Rankine-Hugoniot conditions, from which the name here is borrowed — are the system version of (8), one equation per conserved quantity. They are what determines the pressure, density and temperature behind a blast wave or a bow shock from the conditions ahead of it and the shock speed, without any knowledge of the internal structure of the front. The entropy inequalities are, in that setting, the statement that the thermodynamic entropy of a fluid element increases as it crosses the front, which is where the name comes from; a jump violating them would be a rarefaction shock, never observed.

River and glacier waves, and why the front steepens. Any model in which the propagation speed increases with the amplitude produces the same phenomenon: a flood wave in a channel steepens at its front and decays at its back, and a wave on a glacier or a dune field does the same. The breaking-time formula of this unit gives the time scale on which that steepening is complete, and the Rankine-Hugoniot speed gives the velocity of the resulting front — which, note, is neither the speed of the material behind it nor that in front of it, but a weighted mean determined by the flux.

Numerics. The reason finite-difference schemes for these equations are written in conservative form, as updates of cell averages by numerical fluxes at cell faces, is Step 2 above: a scheme that discretises $u_t + a(u)u_x = 0$ rather than $u_t + F(u)_x = 0$ satisfies no discrete analogue of (9), and will converge to a jump travelling at the wrong speed, even as the mesh is refined. The entropy condition plays the matching role of excluding the spurious stationary jumps such a scheme can otherwise admit. See finite differences for the discretisations and von Neumann stability analysis for the separate question of stability, which is a condition on the time step, not on the form of the scheme.

Reference table

Types of boundary and initial condition

Adapted

Add to PDF

Defines the Dirichlet, Neumann, Robin, mixed, Cauchy and periodic conditions, tabulates which of them give a unique and stable solution for elliptic, parabolic and hyperbolic equations on open and on closed boundaries, and explains each verdict — with a proof of uniqueness for the elliptic case and the compatibility condition the Neumann problem needs.

Prerequisites: Classification of second-order equations, Green's identities, The divergence theorem

Definition(The six kinds of condition)AdaptedAdapted from Strauss (2008), §1.4, pp. 20–22, and §4.1 Exercises 1 and 4, p. 92; Arfken, Weber & Harris (2013), §9.3, pp. 411–412; Hassani (2013), Definitions 21.1.7–21.1.8, pp. 642–643

Let the equation hold in a domain $D$ with boundary $\partial D$ , let $\mathbf{n}$ be the outward unit normal and write $\partial u/\partial n = \mathbf{n}\cdot\nabla u$ for the outward normal derivative. On $\partial D$ , or on the part of it under discussion, the condition is called

Name	Condition	Prescribed data
Dirichlet	$u = g$	the value of $u$
Neumann	$\partial u/\partial n = g$	the normal derivative of $u$
Robin	$\partial u/\partial n + a\,u = g$	a combination, with $a$ a given function
Mixed	Dirichlet on one part of $\partial D$ , Neumann (or Robin) on another	one condition at each boundary point
Cauchy	$u = g$ and $\partial u/\partial n = h$	both, at the same points
Periodic	$u$ and $\partial u/\partial n$ agree at identified boundary points	nothing; the domain is closed up

A condition is homogeneous if the prescribed function vanishes, and inhomogeneous otherwise. An initial condition is a boundary condition imposed at one end of the range of the time variable: $u(x,0) = \varphi(x)$ for the diffusion and Schrödinger equations, and the pair $u(x,0) = \varphi(x)$ , $u_t(x,0) = \psi(x)$ for the wave equation — which is Cauchy data on the surface $t = 0$ . Dirichlet and Neumann conditions are the two cases in which Cauchy data are split; Robin interpolates between them, reducing to Neumann when $a \equiv 0$ and approaching Dirichlet as $a \to \infty$ with $g/a$ held fixed. Hassani defines the Dirichlet and Neumann boundary-value problems specifically as those in which the datum is given on a closed hypersurface, which is the case in which they are well posed for elliptic equations.

Physically, for a vibrating string held at its ends the condition is homogeneous Dirichlet; a free end sliding without resistance gives $u_x = 0$ , a Neumann condition; an end attached to a spring obeying Hooke's law gives a Robin condition, through which the string exchanges energy with the spring; and an end moved in a prescribed way gives an inhomogeneous Dirichlet condition. For heat flow, a perfectly insulated boundary is Neumann with $g = 0$ , a boundary in perfect contact with a reservoir at temperature $g$ is Dirichlet, and Newton cooling into an ambient temperature is Robin. In electrostatics, Dirichlet prescribes the potential on a conductor and Neumann prescribes the normal field, hence the surface charge density. Periodic conditions arise when the domain is physically a ring: diffusion in a closed circular tube of circumference $2l$ , with arclength $-l \le x \le l$ , carries the conditions $u(-l,t) = u(l,t)$ and $u_x(-l,t) = u_x(l,t)$ .

Note(Which conditions suit which equation)AdaptedAdapted from Arfken, Weber & Harris (2013), §9.3, Table 9.1 and the caveats, pp. 412–413; Riley, Hobson & Bence (2006), §20.6.2, Table 20.1 and the uniqueness-to-within-a-constant qualification, p. 705; Hassani (2013), Box 21.1.9, p. 643; Strauss (2008), §4.3, Eqs. (1)–(2) and the sign discussion, pp. 92–93

The following is Arfken, Weber and Harris's table, for the three classes of second-order equation in two variables, with the prototypes Laplace or Poisson in $(x,y)$ , the wave equation in $(x,t)$ and the diffusion equation in $(x,t)$ . "Open" and "closed" refer to the boundary on which the data are given: a curve with endpoints, or a curve enclosing the region. Part of a closed boundary may lie at infinity, provided $u$ or $\partial u/\partial n$ is prescribed there.

Conditions	Boundary	Elliptic	Hyperbolic	Parabolic
Cauchy	open	ill posed (unstable)	well posed	over-determined
Cauchy	closed	over-determined	over-determined	over-determined
Dirichlet	open	under-determined	under-determined	well posed, forwards in time only
Dirichlet	closed	well posed	many solutions	over-determined
Neumann	open	under-determined	under-determined	well posed, forwards in time only
Neumann	closed	well posed	many solutions	over-determined

The entries have been relabelled for compactness. Arfken, Weber and Harris write, in place of "well posed", that the solution is unique and stable; in place of "over-determined", that the conditions are too restrictive; in place of "under-determined", that they are insufficient; in place of "ill posed (unstable)", that the results are unphysical because of instability; in place of "many solutions", that the solution is not unique; and in place of "forwards in time only", that the solution is unique and stable in one direction.

Reading off only the well-posed rows leaves three pairings, which are Riley, Hobson and Bence's Table 20.1 and, in $m$ variables, Hassani's Box 21.1.9. A hyperbolic equation takes Cauchy data on an open boundary. A parabolic equation takes one condition of Dirichlet or Neumann type, again on an open boundary. An elliptic equation takes one condition of Dirichlet or Neumann type on a closed boundary.

Robin and mixed conditions impose, like Dirichlet and Neumann, one scalar condition at each boundary point. For a mixed condition — Dirichlet on part of the boundary, Neumann or Robin on the rest — the table's verdicts carry over unchanged. For Robin the sign of the coefficient matters, and Strauss is explicit about why. He writes the one-dimensional conditions as $X' - a_0X = 0$ at $x = 0$ and $X' + a_lX = 0$ at $x = l$ , with the two constants carrying opposite signs because the outward normal points left at one end and right at the other, and he reads $a_0, a_l > 0$ as radiation of energy, $a_0, a_l < 0$ as absorption and $a_0 = a_l = 0$ as insulation. He warns at the outset that, because of that asymmetry in the normal, each of the two constants may affect the eigenfunctions in the opposite way to the other. In the radiating case — $a \ge 0$ in the convention $\partial u/\partial n + a u = g$ — the elliptic verdicts do carry over, because the boundary term in the uniqueness proof below then has a sign. A sufficiently negative coefficient does not: it can put an eigenvalue of the associated problem at zero, and uniqueness on a closed boundary then fails. Example (a) below exhibits such a case, and shows what $a \ge 0$ buys. Periodic conditions are the case in which the boundary has been removed altogether by identifying its two parts, so for a parabolic or hyperbolic equation on a ring only the initial conditions remain to be given.

Both sources warn that the table simplifies. Riley, Hobson and Bence say the general question of which conditions are appropriate is complex and do not pursue it, and they qualify "appropriate" as giving a solution that is unique — sometimes only to within an additive constant — and well defined. Arfken, Weber and Harris call the part played by boundary conditions too complicated to analyse exactly within the scope of their book, say that their Table 9.1 oversimplifies in several respects, and give the Helmholtz equation $\nabla^2\psi \pm k^2\psi = 0$ as an example the table mis-describes: with Dirichlet conditions on a closed boundary it has solutions only for certain values of $k$ , those being the eigenvalues (see the Helmholtz equation).

Remark(Why each verdict holds)Original

The entries are not independent; most follow from two facts about characteristics and one counting principle.

Elliptic, closed, Dirichlet or Neumann: unique and stable. The uniqueness theorem below says so, and its proof turns the problem into $\int_D|\nabla w|^2 = 0$ for the difference $w$ of two solutions. The elliptic Dirichlet and Neumann problems are therefore the prototype well-posed problems, as Stakgold and Holst put it for $-\Delta u = q$ in $\Omega$ with $u = f$ on the boundary.
Elliptic, closed, Cauchy: too restrictive. Dirichlet data alone already fix the solution; the normal derivative on the boundary is then determined, not free. Prescribing it as well generally contradicts what the Dirichlet data imply. Arfken, Weber and Harris make exactly this point: since Dirichlet and Neumann conditions each give a unique solution, and the two solutions are different, Cauchy data — which are both at once — invite an inconsistency.
Elliptic, open, Cauchy: unstable. This is Hadamard's example, worked out in the classification unit: for Laplace's equation in a half-plane there are solutions whose Cauchy data on the boundary line tend to zero while the solutions themselves blow up at every interior point. The absence of real characteristics, which one might expect to make the Cauchy problem harmless, does not save it.
Elliptic, open, Dirichlet or Neumann: insufficient. An open arc leaves the solution unconstrained on the rest of the boundary, and an elliptic equation propagates nothing along curves; so the data simply do not reach the interior.
Hyperbolic, open, Cauchy: unique and stable. Here there are two families of real characteristics. Cauchy data on an arc $C$ that is nowhere characteristic are carried along the characteristic of each family through each point of $C$ , and the solution is determined in the region covered by both families. For the wave equation $u_{xx} - c^{-2}u_{tt} = 0$ the characteristics are $x \pm ct = \text{const}$ , and data on the segment $t = 0$ , $0 \le x \le L$ determine $u$ exactly in the triangle those lines cut out (see d'Alembert's solution).
Hyperbolic, closed, Dirichlet or Neumann: not unique; closed, Cauchy: too restrictive. One datum per point is the wrong count for an equation needing two on an initial surface, and a closed boundary generally meets the same characteristic twice. Riley, Hobson and Bence note that a boundary curve crossing a characteristic more than once can over-determine the problem so that no solution exists, and that the cure is to relax the conditions on the offending parts of the curve from Cauchy to Dirichlet or Neumann.
Parabolic, open, Dirichlet or Neumann: unique and stable in one direction. The standard parabolic problem is Evans's: the equation in $U\times(0,T]$ , one condition on the lateral boundary $\partial U\times[0,T]$ , and one initial condition on $U\times\{t=0\}$ . In the $x$ – $t$ plane this data set is an open curve — the bottom of a rectangle and its two sides — and it determines $u$ on the remaining, open side. "In one direction" means forward in time only: the backward problem is ill posed, because diffusion smooths, so reversing it amplifies small, high-frequency errors without bound.
Parabolic, anything, Cauchy: too restrictive; Dirichlet or Neumann on a closed boundary: too restrictive. A parabolic equation is first order in time, so one condition per boundary point and one initial function are already enough; adding a terminal condition at $t = T$ , which closing the boundary does, over-determines it.

Reference table(Uniqueness for the elliptic Dirichlet and Neumann problems)AdaptedAdapted from Riley, Hobson & Bence (2006), §20.7, Uniqueness theorem, p. 706

Let $V$ be a connected region with boundary $S$ . (Riley, Hobson and Bence state the theorem for "a region $V$ and its boundary $S$ " and do not say connected; the word is added here because the last step of their proof — $\nabla w = 0$ , hence $w$ constant — needs it, as the remark below records.) Suppose $u$ is real-valued, that its first and second partial derivatives all exist and are continuous on $V \cup S$ , and that

\nabla^2 u = \rho \ \text{ in } V, \qquad\text{together with either}\quad u = f \ \text{ on } S \quad\text{or}\quad \frac{\partial u}{\partial n} = g \ \text{ on } S ,

where $\rho$ , $f$ and $g$ are prescribed. Then two functions satisfying all of this can differ by at most an additive constant. Under the Dirichlet condition that constant is forced to be zero, so the solution is unique; under the Neumann condition it is free, and the solution is unique only modulo constants, which is how Riley, Hobson and Bence state the theorem.

Idea of the proofOriginal

Subtracting two solutions kills both the source and the data, leaving a harmonic function whose boundary value, or whose normal derivative, vanishes. Green's first identity converts the volume integral of $|\nabla w|^2$ into a boundary integral of $w\,\partial w/\partial n$ plus a volume integral of $w\nabla^2 w$ . Each of the two boundary conditions kills the boundary term, and harmonicity kills the other; so the integral of a non-negative quantity vanishes, and the gradient must vanish identically.

ProofAdaptedAdapted from Riley, Hobson & Bence (2006), §20.7, p. 706

Suppose $u_1$ and $u_2$ both satisfy the hypotheses, and set $w = u_1 - u_2$ .

Step 1 (the difference is harmonic). Subtracting the two equations, $\nabla^2 w = \rho - \rho = 0$ in $V$ .

Step 2 (the difference has vanishing data). In the Dirichlet case $u_1 = f = u_2$ on $S$ , so $w = 0$ there; in the Neumann case $\partial u_1/\partial n = g = \partial u_2/\partial n$ , so $\partial w/\partial n = 0$ there. The two cases are not mixed: whichever datum is prescribed, the corresponding quantity vanishes for $w$ .

Step 3 (Green's first identity). Taking both scalar fields in Green's first identity to be $w$ , which the assumed smoothness of $u_1$ and $u_2$ on $V\cup S$ permits,

\int_V \Bigl( w\,\nabla^2 w + |\nabla w|^2 \Bigr)\,dV = \int_S w\,\frac{\partial w}{\partial n}\,dS .

Step 4 (both unwanted terms vanish). The first term on the left vanishes by Step 1. The right-hand side vanishes by Step 2: the integrand is $w\,\partial w/\partial n$ , and one of its two factors is identically zero on $S$ in either case. Hence

\int_V |\nabla w|^2\,dV = 0 .

Step 5 (conclude). The integrand is continuous and non-negative, so it vanishes identically: $\nabla w = 0$ throughout $V$ . Hence $w$ takes a single value on each connected component of $V$ , so $u_1$ and $u_2$ differ by a constant. With Dirichlet data $w$ is zero on $S$ , which pins the constant at zero and gives $u_1 \equiv u_2$ . With Neumann data nothing pins it, and two solutions may differ by any constant whatever. $\square$

RemarkOriginal

Where the hypotheses are used, and the Neumann problem's extra condition.

Smoothness up to the boundary is used in Step 3; Green's first identity needs $w$ differentiable in $V$ and integrable on $S$ . Without it the boundary integral need not exist and the argument has no first line.
Realness is used in Step 5, where $|\nabla w|^2 \ge 0$ is the whole argument. For a complex-valued $w$ one applies the same identity to $\bar w$ and $w$ instead.
Connectedness is used in Step 5. On a region with two components, $w$ may be a different constant on each; with Dirichlet data each component touches $S$ , so each constant is zero.
Uniqueness is not existence. The theorem says a solution, if found by any means at all, is the only one; this is what licenses the method of images and every other trick that guesses a solution and verifies it (see the method of images).
The Neumann problem needs a compatibility condition. For a solution of $\nabla^2 u = \rho$ in $V$ with $\partial u/\partial n = g$ on $S$ to exist at all, integrate the equation over $V$ and apply the divergence theorem to $\nabla u$ :
$\int_V \rho \,dV = \int_V \nabla\cdot(\nabla u)\,dV = \int_S \frac{\partial u}{\partial n}\,dS = \int_S g \,dS .$
So $\int_V\rho\,dV = \int_S g\,dS$ is necessary. If the data violate it, there is no solution whatever, and the table's verdict "unique, stable solution" for the elliptic Neumann problem must be read with that proviso and with the arbitrary constant of Step 5. Strauss sets precisely this calculation, for $\nabla^2 u = f$ with $\partial u/\partial n = 0$ , as an exercise, and asks the reader for a physical interpretation without giving one; the following reading is therefore ours. For steady heat flow with $\rho$ a source density, the displayed identity says that the total heat generated inside $V$ must equal the total flux out through $S$ , as it must if nothing is to accumulate. In Strauss's own case the boundary is insulated, so the flux out is zero and the identity reduces to the statement that the net heat generated must vanish.

ExampleOriginal

(a) A Robin problem, solved. Steady heat conduction in a slab $0 < x < 1$ with no internal source satisfies $u'' = 0$ . Suppose both faces exchange heat with surroundings, giving the Robin conditions

u'(0) - 2u(0) = -2, \qquad u'(1) + u(1) = 4 .

Every solution of $u''=0$ is $u = A + Bx$ , so $u' = B$ and the conditions read $B - 2A = -2$ and $B + (A + B) = 4$ , that is $A + 2B = 4$ . From the first, $B = 2A - 2$ ; substituting, $A + 4A - 4 = 4$ , so $A = \tfrac85$ and $B = \tfrac65$ :

u(x) = \frac{8}{5} + \frac{6}{5}x , \qquad u(0) = 1.6, \quad u(1) = 2.8 .

Checking: $u'(0) - 2u(0) = 1.2 - 3.2 = -2$ and $u'(1) + u(1) = 1.2 + 2.8 = 4$ . Two conditions, one at each end of an interval whose boundary is the pair of endpoints — a closed boundary in one dimension — give a unique solution, as the elliptic row of the table predicts.

Why the sign of the Robin coefficient was needed. In the convention $\partial u/\partial n + a u = g$ the outward normal is $-1$ at $x = 0$ and $+1$ at $x = 1$ , so the two conditions above read $\partial u/\partial n + 2u = 2$ at the left end and $\partial u/\partial n + u = 4$ at the right: both have $a > 0$ , Strauss's radiating case. Take instead $a = -2$ at both ends with zero data, that is $u'' = 0$ on $(0,1)$ with

-u'(0) - 2u(0) = 0, \qquad u'(1) - 2u(1) = 0 .

Then $u(x) = 1 - 2x$ satisfies both: $u' \equiv -2$ , so the left condition is $2 - 2\cdot 1 = 0$ and the right is $-2 - 2\cdot(-1) = 0$ . Since $u \equiv 0$ also satisfies them, uniqueness fails, and so the inhomogeneous problem is not solvable for every $g$ . In Strauss's notation this is $a_0 = a_l = -2$ , his absorbing case, and $1 - 2x$ is an eigenfunction of $-X'' = \lambda X$ with those conditions at eigenvalue $\lambda = 0$ . With $a \ge 0$ nothing of the sort can happen. Running Steps 1 to 4 of the proof on the difference $w$ of two solutions, with $\partial w/\partial n = -a w$ on $S$ , gives

\int_V |\nabla w|^2\,dV = \int_S w\,\frac{\partial w}{\partial n}\,dS = -\int_S a\,w^2\,dS ,

whose left side is non-negative and whose right side is non-positive; both therefore vanish, so $w$ is constant, and the constant is zero wherever $a > 0$ on $S$ .

(b) The Neumann obstruction, in one dimension. Replace the Robin conditions by $u'(0) = u'(1) = 0$ and allow a source: $u'' = f$ on $(0,1)$ . Integrating once gives $u'(1) - u'(0) = \int_0^1 f\,dx$ , so the compatibility condition is $\int_0^1 f \,dx = 0$ — exactly the identity of the remark, with $S$ the two endpoints.

Incompatible data. For $f \equiv 1$ we would need $0 = 1$ . No solution exists, however smooth everything is.
Compatible data. For $f(x) = \cos 2\pi x$ we have $\int_0^1 f = 0$ . Integrating, $u' = \sin(2\pi x)/(2\pi)$ , which vanishes at both ends as required, and $u(x) = -\cos(2\pi x)/(4\pi^2) + D$ for every constant $D$ . So there are solutions, and they form a one-parameter family: uniqueness holds only up to the additive constant, as Step 5 says. Here $u(0) = D - 1/(4\pi^2) = D - 0.025330$ and $u(\tfrac12) = D + 0.025330$ .

(c) Counting conditions on the three prototypes. On $0 < x < L$ :

Equation	Conditions needed	A standard choice
$u_{xx} + u_{yy} = 0$ in a rectangle	one per boundary point, all round	$u$ given on all four sides
$u_t = k u_{xx}$	one initial function, one condition at each end	$u(x,0) = \varphi$ , $u(0,t) = u(L,t) = 0$
$u_{tt} = c^2u_{xx}$	two initial functions, one condition at each end	$u(x,0) = \varphi$ , $u_t(x,0) = \psi$ , $u(0,t) = u(L,t) = 0$

The wave equation needs the extra initial function because it is second order in $t$ ; the mixed choice $u(0,t) = u_x(L,t) = 0$ for either evolution equation is legitimate and, for the diffusion equation, produces the eigenfunctions $\sin\bigl[(n+\tfrac12)\pi x/L\bigr]$ , $n = 0, 1, 2, \dots$ , instead of the pure sines $\sin(n\pi x/L)$ of the Dirichlet problem. These do satisfy both conditions: the sine vanishes at $x = 0$ , and its derivative at $x = L$ is proportional to $\cos\bigl[(n+\tfrac12)\pi\bigr] = 0$ . The eigenvalues are $\bigl[(n+\tfrac12)\pi/L\bigr]^2$ : the half-integer multiples of $\pi/L$ , squared, in place of the integer ones. Swapping the two ends, to $u_x(0,t) = u(L,t) = 0$ , replaces these by the cosines $\cos\bigl[(n+\tfrac12)\pi x/L\bigr]$ with the same eigenvalues. On a ring of circumference $2l$ the periodic conditions $u(-l,t) = u(l,t)$ , $u_x(-l,t) = u_x(l,t)$ replace the two end conditions, and the eigenvalues become $(n\pi/l)^2$ for $n = 0, 1, 2, \dots$ , with the value $n = 0$ (a constant eigenfunction) present and every $n \ge 1$ doubled, since both $\cos(n\pi x/l)$ and $\sin(n\pi x/l)$ qualify.

ApplicationOriginal

Choosing the condition is the modelling step that decides whether a computation means anything, and the table is the checklist.

In a thermal or electrostatic finite-element calculation, Dirichlet conditions represent contact with a reservoir or a held potential, Neumann conditions represent insulation or a prescribed flux, and Robin conditions represent convective exchange with an ambient medium at a given temperature, with $a$ the heat-transfer coefficient divided by the conductivity; a pure Neumann problem will not converge unless the fluxes balance to within round-off, which is the compatibility condition of the remark appearing as a solver diagnostic, and the solution it returns is pinned only after one node is anchored, which is the additive constant.

In ocean and atmosphere models the lateral walls of a regional domain are open boundaries for a hyperbolic system, so the data must be posed in the characteristic variables that enter the domain and nothing may be imposed on those that leave; imposing a full set of values on all variables there is the closed-boundary Cauchy case of the table, and it produces spurious reflections, the numerical form of "too restrictive". Periodic conditions are used for the zonal direction of a global model, where the domain really is a circle, and in molecular-dynamics and spectral computations to remove boundaries entirely.

For the reverse problems — reconstructing a past temperature field, or an interior source from measurements on an enclosing surface — the table says what is being attempted: the backward parabolic problem and the elliptic Cauchy problem, both ill posed. They are solved only after regularisation, that is, after being replaced by a nearby well-posed problem, which is why such inversions always carry a smoothing parameter.

Method

Similarity and self-similar solutions

Adapted

Add to PDF

Shows how a dilation symmetry of a partial differential equation forces solutions of the form t^{-alpha} v(x/t^beta), reduces the equation to an ordinary differential equation for the profile v, and carries this out for the diffusion equation, for the porous medium equation (the Barenblatt solution) and for the Airy equation.

Prerequisites: Classification of second-order equations, The heat kernel, Partial derivatives and the chain rule

Definition(Scaling symmetry, similarity variable, self-similar solution)AdaptedAdapted from Evans (2010), §2.3.1a, p. 43, and §4.2.2, pp. 175–176; Strauss (2008), §2.4, property (e), p. 47

Let an equation be posed for $u(x,t)$ with $x \in \mathbb{R}^n$ , $t > 0$ . Say that the equation has a dilation symmetry with exponents $(\alpha, \beta)$ if, whenever $u$ is a solution, so is

u_\lambda(x,t) := \lambda^{\alpha}\,u\bigl(\lambda^{\beta}x,\ \lambda t\bigr) \qquad\text{for every } \lambda > 0 . \tag{1}

A solution is self-similar, or a similarity solution, if it is unchanged by the whole family: $u = u_\lambda$ for all $\lambda > 0$ . Setting $\lambda = t^{-1}$ in that identity and writing $v(y) := u(y,1)$ gives the general form of such a solution,

u(x,t) = \frac{1}{t^{\alpha}}\,v\!\left(\frac{x}{t^{\beta}}\right) , \tag{2}

in which

y = \frac{x}{t^{\beta}}

is the similarity variable and $v$ the profile. The profile is a function of $n$ variables instead of $n+1$ : the symmetry has been traded for a reduction in dimension. When $v$ is additionally radial, $v(y) = w(|y|)$ , the reduction is to a single ordinary differential equation.

The name says what the picture is. At two different times the graph of $u$ is the same shape, rescaled: widened by the factor $(t_2/t_1)^{\beta}$ in $x$ and shortened by $(t_1/t_2)^{\alpha}$ in amplitude. Strauss states the symmetry of the diffusion equation $u_t = ku_{xx}$ in the equivalent form that $u(\sqrt a\,x, at)$ solves it whenever $u$ does, for every $a > 0$ , which is (1) with $\alpha = 0$ , $\beta = \tfrac12$ and $\lambda = a$ ; the chain rule gives it at once, since the dilation multiplies $u_t$ by $a$ and $u_{xx}$ by $(\sqrt a)^2 = a$ .

Method(Finding a similarity solution)AdaptedAdapted from Evans (2010), §2.3.1a, pp. 43–44, and §4.2.2, pp. 175–177; Stakgold & Holst (2011), Exercise 8.2.6(c), p. 488

Find the dilation symmetry. Substitute $u_\lambda$ from (1) into the equation and collect the powers of $\lambda$ that multiply each term. The equation is invariant exactly when all those powers agree; this gives one linear relation between $\alpha$ and $\beta$ . (Equivalently: count the derivatives. One $t$ -derivative scales like $\lambda$ , one $x$ -derivative like $\lambda^{\beta}$ , so $u_t$ and $\Delta u$ balance when $\beta = \tfrac12$ , which is the diffusive scaling $x \sim \sqrt{t}$ .)
Fix the remaining exponent from the side conditions. One relation leaves a one-parameter family of candidate scalings. Choose within it by asking that the data be invariant too. Three cases recur:
- a conserved quantity, usually $\int u \,dx$ , forces $\alpha = n\beta$ ;
- initial or boundary data that are themselves dilation-invariant — a step, a constant, a delta — force the exponent that leaves them unchanged;
- a radial reduction that is to simplify the resulting ordinary differential equation may itself dictate the choice.
Write the ansatz and substitute. Put $u = t^{-\alpha}v(x\,t^{-\beta})$ into the equation and use the chain rule. With $y = x t^{-\beta}$ the time derivative contributes $-t^{-(\alpha+1)}\bigl(\alpha v + \beta\,y\cdot Dv\bigr)$ , and each pair of space derivatives contributes a factor $t^{-2\beta}$ .
Check that the powers of $t$ cancel. If steps 1 and 2 were done correctly, every term carries the same power of $t$ ; divide it out. What is left is an equation for $v$ in the variable $y$ alone.
Reduce to an ordinary differential equation and integrate it. If the equation and the data are rotation-invariant, look for $v(y) = w(r)$ with $r = |y|$ ; the reduced equation is then an ordinary differential equation in $r$ , usually of the form $(r^{n-1}\,\cdot)' = 0$ after one regrouping, so it integrates once immediately. Decay as $r \to \infty$ sets the first constant of integration to zero.
Fix the remaining constants and the domain. Impose the normalisation (total mass, or a boundary value) and check the sign of $w$ ; where the formula would go negative for a non-negative quantity, the solution is the positive part, and the boundary of its support is a free boundary.
Verify. Substitute the closed form back into the original equation. This step is not optional: steps 3 and 4 are easy to get wrong by one power.

Example(The diffusion similarity solution)AdaptedAdapted from Strauss (2008), §2.4, Steps 1–4, pp. 47–49

Solve $u_t = k u_{xx}$ on $-\infty < x < \infty$ , $t > 0$ , with the step initial datum $u(x,0) = 1$ for $x > 0$ and $0$ for $x < 0$ . Call the solution $Q$ .

Symmetry, and why the data pick it out. The equation admits $u(x,t) \mapsto u(\sqrt a\,x, at)$ for every $a > 0$ . The step datum is also unchanged by $x \mapsto \sqrt a\,x$ , because only the sign of $x$ matters. So the problem as a whole is dilation-invariant, and the solution it determines must be too. A function of $(x,t)$ invariant under $(x,t)\mapsto(\sqrt a\,x, at)$ for all $a>0$ depends on $x$ and $t$ only through the combination $x/\sqrt{t}$ , since that combination is the invariant one. This is step 2 of the method in the form "invariant data".
Ansatz. Write $Q(x,t) = g(p)$ with $p = x/\sqrt{4kt}$ , the constant $4k$ inserted only to tidy the algebra.
Reduce. By the chain rule, $\partial p/\partial t = -p/(2t)$ and $\partial p/\partial x = 1/\sqrt{4kt}$ , so
$Q_t = -\frac{p}{2t}\,g'(p), \qquad Q_x = \frac{g'(p)}{\sqrt{4kt}}, \qquad Q_{xx} = \frac{g''(p)}{4kt} .$
Hence $0 = Q_t - kQ_{xx} = -\frac{1}{t}\Bigl(\frac{p}{2}g' + \frac14 g''\Bigr)$ , and multiplying by $-4t$ ,
$g'' + 2p\,g' = 0 .$
Every power of $t$ has cancelled, as step 4 requires.
Integrate. The integrating factor is $e^{p^2}$ , since $(e^{p^2}g')' = e^{p^2}(g'' + 2pg')$ . So $g'(p) = c_1 e^{-p^2}$ and
$g(p) = c_1\int_0^p e^{-q^2}\,dq + c_2 .$
Match the data. Letting $t \downarrow 0$ with $x > 0$ sends $p \to +\infty$ , and with $x < 0$ sends $p \to -\infty$ . Using $\int_0^{\pm\infty}e^{-q^2}dq = \pm\sqrt\pi/2$ , the two requirements $1 = c_1\sqrt\pi/2 + c_2$ and $0 = -c_1\sqrt\pi/2 + c_2$ give $c_1 = 1/\sqrt\pi$ and $c_2 = \tfrac12$ . Therefore
$Q(x,t) = \frac12 + \frac{1}{\sqrt\pi}\int_0^{x/\sqrt{4kt}} e^{-q^2}\,dq = \frac12\left[1 + \operatorname{erf}\!\left(\frac{x}{\sqrt{4kt}}\right)\right] \qquad (t > 0) .$
Differentiate to get the kernel. Since the equation is linear with constant coefficients, $\partial Q/\partial x$ is again a solution, and
$S(x,t) = \frac{\partial Q}{\partial x} = \frac{1}{\sqrt{4\pi k t}}\,e^{-x^2/4kt} ,$
which is the heat kernel; convolving it with any initial datum $\varphi$ solves the general initial-value problem. Note that $S$ is itself self-similar, with $\alpha = \tfrac12$ and $\beta = \tfrac12$ in the notation of (2): $S(x,t) = t^{-1/2}\,\sigma(x/\sqrt t)$ with $\sigma(y) = (4\pi k)^{-1/2}e^{-y^2/4k}$ , and $\alpha = n\beta$ with $n = 1$ , which is the conserved-mass condition of step 2.
Numbers. For $k = 1$ : $Q(1,1) = \tfrac12\bigl[1 + \operatorname{erf}(\tfrac12)\bigr] = \tfrac12(1 + 0.5204999) = 0.7602499$ , and $S(1,1) = (4\pi)^{-1/2}e^{-1/4} = 0.2196956$ . Substituting the closed form for $Q$ into $Q_t - kQ_{xx}$ gives zero identically.
The same method on a half-line. Stakgold and Holst treat $w_t = w_{xx}$ for $x > 0$ with $w(x,0) = 0$ and $w(0,t) = 1$ , which is invariant under $x^* = ax$ , $t^* = a^2t$ ; taking $a = 1/(2\sqrt t)$ turns the invariance into $w(x,t) = h\bigl(x/2\sqrt t\bigr)$ , the reduced equation is $h'' + 2\xi h' = 0$ with $h(0) = 1$ , $h(\infty) = 0$ , and the solution is $h = 1 - \operatorname{erf}\xi$ , so $w(x,t) = 1 - \operatorname{erf}\bigl(x/2\sqrt t\bigr)$ . The same reduced equation appears, with different boundary conditions picking a different pair of constants.

Example(The Barenblatt solution of the porous medium equation)AdaptedAdapted from Evans (2010), §4.2.2, Eqs. (27)–(36), pp. 175–177

The porous medium equation is

u_t - \Delta\bigl(u^{\gamma}\bigr) = 0 \quad\text{in } \mathbb{R}^n\times(0,\infty), \qquad u \ge 0,\ \gamma > 1 ,

a nonlinear diffusion in which the diffusivity grows with $u$ . It is quasilinear and degenerate: where $u = 0$ the diffusion switches off.

Symmetry. Substituting $u_\lambda = \lambda^{\alpha}u(\lambda^{\beta}x, \lambda t)$ , the term $u_t$ acquires $\lambda^{\alpha+1}$ while $\Delta(u^{\gamma})$ acquires $\lambda^{\alpha\gamma + 2\beta}$ . The two balance exactly when
$\alpha + 1 = \alpha\gamma + 2\beta . \tag{3}$
Reduce. With $u = t^{-\alpha}v(y)$ , $y = t^{-\beta}x$ , the chain rule gives
$\alpha\,t^{-(\alpha+1)}v + \beta\,t^{-(\alpha+1)}\,y\cdot Dv + t^{-(\alpha\gamma+2\beta)}\Delta\bigl(v^{\gamma}\bigr) = 0 ,$
so (3) makes every power of $t$ equal and leaves $\alpha v + \beta\,y\cdot Dv + \Delta(v^{\gamma}) = 0$ .
Go radial. Put $v(y) = w(r)$ , $r = |y|$ . Using $\Delta f(r) = f'' + \frac{n-1}{r}f'$ ,
$\alpha w + \beta r w' + \bigl(w^{\gamma}\bigr)'' + \frac{n-1}{r}\bigl(w^{\gamma}\bigr)' = 0 .$
Choose the second exponent to make it integrable. Multiply by $r^{n-1}$ . The last two terms become $\bigl(r^{n-1}(w^{\gamma})'\bigr)'$ , and the first two become $\alpha r^{n-1}w + \beta r^n w'$ , which is the derivative $\beta\bigl(r^n w\bigr)' = \beta\bigl(n r^{n-1}w + r^n w'\bigr)$ precisely when
$\alpha = n\beta . \tag{4}$
This is the conserved-mass condition of step 2 of the method, and here it is also exactly what makes the equation a total derivative:
$\bigl(r^{n-1}(w^{\gamma})'\bigr)' + \beta\bigl(r^n w\bigr)' = 0 .$
Integrate twice. Hence $r^{n-1}(w^{\gamma})' + \beta r^n w = a$ for a constant $a$ , and $a = 0$ if $w, w' \to 0$ fast enough as $r \to \infty$ . Dividing by $r^{n-1}$ leaves $(w^{\gamma})' = -\beta r w$ . Writing $(w^{\gamma})' = \frac{\gamma}{\gamma-1}\,w\,(w^{\gamma-1})'$ , which holds wherever $w > 0$ , this becomes $(w^{\gamma-1})' = -\frac{\gamma-1}{\gamma}\beta r$ , so
$w^{\gamma-1} = b - \frac{\gamma-1}{2\gamma}\beta r^2$
for a constant $b$ , and taking the positive part to keep $w \ge 0$ ,
$w(r) = \left(b - \frac{\gamma-1}{2\gamma}\beta r^2\right)_{+}^{1/(\gamma-1)} .$
Solve (3) and (4) for the exponents. Substituting $\alpha = n\beta$ into (3) gives $n\beta + 1 = n\beta\gamma + 2\beta$ , so $\beta\bigl(n(\gamma-1) + 2\bigr) = 1$ and
$\alpha = \frac{n}{n(\gamma-1)+2}, \qquad \beta = \frac{1}{n(\gamma-1)+2} .$
Result. Undoing the substitutions,
$u(x,t) = \frac{1}{t^{\alpha}}\left(b - \frac{\gamma-1}{2\gamma}\,\beta\,\frac{|x|^2}{t^{2\beta}}\right)_{+}^{1/(\gamma-1)} \qquad (x \in \mathbb{R}^n,\ t>0) ,$
the Barenblatt-Kompaneetz-Zeldovich solution, or Barenblatt solution, with $b > 0$ a free parameter fixing the total mass.
Check, with numbers. For $n = 1$ , $\gamma = 2$ the exponents are $\alpha = \beta = \tfrac13$ and the formula reads $u = t^{-1/3}\bigl(b - x^2/(12\,t^{2/3})\bigr)_+$ . Substituting into $u_t = (u^2)_{xx}$ inside the support confirms it, as does the case $n = 2$ , $\gamma = 3$ , where $\alpha = \tfrac13$ , $\beta = \tfrac16$ . With $b = 1$ , $n = 1$ , $\gamma = 2$ : at $t = 1$ the support is $|x| \le \sqrt{12} = 3.4641$ and the peak is $u(0,1) = 1$ ; at $t = 8$ the support has widened to $|x| \le \sqrt{12}\cdot 8^{1/3} = 6.9282$ and the peak has dropped to $8^{-1/3} = 0.5$ . The total mass $\int u\,dx$ is $\tfrac43\sqrt{12}\,b^{3/2}$ at both times, as $\alpha = n\beta$ guarantees.
Why this matters. The Barenblatt solution has compact support at each time: the region $\{u > 0\}$ is a ball whose radius grows like $t^{\beta}$ , with $\beta < \tfrac12$ for $\gamma > 1$ . Evans records that the same holds generally: start a non-negative weak solution of this equation from initial data supported in a bounded set, and the set $\{u>0\}$ stays bounded and spreads at finite speed — in contrast to the linear heat equation, whose kernel is positive everywhere for every $t>0$ , so that a disturbance is felt instantly at arbitrary distance. For that reason the porous medium equation is sometimes preferred as a model of diffusive spreading.

ExampleOriginal

A dispersive similarity solution: the Airy equation from $u_t + u_{xxx} = 0$ . The method is not confined to diffusion. Take the linearised Korteweg-de Vries, or Airy, equation on the line, $u_t + u_{xxx} = 0$ , and look for the self-similar solution of unit mass.

Symmetry. Under $u_\lambda = \lambda^{\alpha}u(\lambda^{\beta}x, \lambda t)$ the term $u_t$ picks up $\lambda^{\alpha+1}$ and $u_{xxx}$ picks up $\lambda^{\alpha+3\beta}$ , so invariance requires $\beta = \tfrac13$ : the natural scaling here is $x \sim t^{1/3}$ , not $x\sim t^{1/2}$ .
Second exponent. The equation conserves mass, because $\frac{d}{dt}\int u\,dx = -\int u_{xxx}\,dx = 0$ for $u$ decaying at infinity. With $n = 1$ , step 2 of the method gives $\alpha = n\beta = \tfrac13$ .
Reduce. Write $u = t^{-1/3}v(\xi)$ with $\xi = x\,t^{-1/3}$ . Then
$u_t = -\tfrac13 t^{-4/3}\bigl(v + \xi v'\bigr), \qquad u_{xxx} = t^{-1/3}\cdot t^{-1}\,v''' = t^{-4/3}v''' ,$
so after dividing by $t^{-4/3}$ the equation becomes
$v''' = \tfrac13\bigl(\xi v\bigr)' .$
Integrate once. Hence $v'' = \tfrac13\xi v + C$ , and $C = 0$ if $v, v''\to 0$ as $\xi \to +\infty$ . So $v'' = \tfrac13 \xi v$ .
Recognise it. Put $s = \xi/3^{1/3}$ , so $d^2/d\xi^2 = 3^{-2/3}d^2/ds^2$ and the equation becomes $3^{-2/3}v_{ss} = 3^{1/3}sv/3 = 3^{-2/3}sv$ , that is
$\frac{d^2v}{ds^2} = s\,v ,$
which is Airy's equation (see Airy functions). The solution decaying as $s \to +\infty$ is $v = c\,\mathrm{Ai}(s)$ ; the second solution $\mathrm{Bi}$ grows and is rejected.
Normalise. With $\int_0^\infty\mathrm{Ai} = \tfrac13$ and $\int_{-\infty}^0 \mathrm{Ai} = \tfrac23$ — the latter an improper integral, convergent only conditionally, since $\mathrm{Ai}(s)$ oscillates with envelope $|s|^{-1/4}$ as $s\to-\infty$ — the total is $\int_{-\infty}^{\infty}\mathrm{Ai}(s)\,ds = 1$ . Since $\int v\,d\xi = c\,3^{1/3}\int \mathrm{Ai}(s)\,ds$ , unit mass needs $c = 3^{-1/3}$ , and
$u(x,t) = \frac{1}{(3t)^{1/3}}\,\mathrm{Ai}\!\left(\frac{x}{(3t)^{1/3}}\right) .$
Check. Differentiating this expression numerically at $x = 0.7$ , $t = 0.4$ gives $u_t = -0.04886036050$ and $u_{xxx} = +0.04886036050$ , which sum to zero to the working precision. At $x = 0$ , $t = \tfrac13$ the value is $\mathrm{Ai}(0) = 3^{-2/3}/\Gamma(2/3) = 0.3550280539$ .
Contrast with diffusion. The profile is oscillatory on one side and exponentially small on the other, and it spreads like $t^{1/3}$ rather than $t^{1/2}$ , with amplitude falling like $t^{-1/3}$ . Nothing in the method changed; only the count of derivatives in step 1 did, and that is what sets the spreading rate.

RemarkOriginal

What the method needs, and what it gives.

A genuine symmetry. Step 1 fails whenever the equation contains a dimensional constant that the dilation does not respect. The equation $u_t - u_{xx} = u$ has no dilation symmetry of the form (1), because the right-hand side scales like $\lambda^{\alpha}$ while $u_t$ scales like $\lambda^{\alpha+1}$ ; the fix is to remove the offending term first, here by substituting $u = e^{t}\tilde u$ . A term with a fixed length or time scale in it always has to be absorbed or scaled away before step 1 can succeed.
Invariant data. The construction produces a solution of the equation, not of an arbitrary initial-value problem. It solves exactly the problems whose data are themselves invariant under the same scaling: a step, a constant boundary value, a point mass. That is why the heat-kernel and Barenblatt cases come out as fundamental solutions, with a delta as the initial condition.
Two conditions, two exponents. One relation comes from the equation (step 1), the second from the data or a conservation law (step 2). If both are supplied by the equation, as happens for some systems, there is no free exponent and the method either works outright or not at all; if neither is, a one-parameter family of similarity solutions results, which is the normal situation for second-kind self-similarity, where the exponent must be found by solving an eigenvalue problem for the profile.
Why these special solutions matter more than they look. Self-similar solutions are usually the long-time attractors of their equations: general initial data forget their shape and approach the similarity profile, with the exponents $\alpha,\beta$ setting the rates. The heat kernel plays that role for the linear diffusion equation, and the Barenblatt profile for the porous medium equation. This is also why dimensional analysis so often predicts the right power law: it is step 1 of this method, done without writing the equation down.
Degeneracy and free boundaries. In step 6 the positive part is not a cosmetic truncation: it is forced, because the equation degenerates where $u = 0$ , and the edge of the support is a genuine moving interface. Taking the positive part of a formula is legitimate only for an equation that admits weak solutions of that kind; doing it to a solution of the linear heat equation would produce something that is not a solution at all.

ApplicationOriginal

Similarity solutions are the standard way of extracting a quantitative law from a diffusive or dispersive process when the geometry is simple enough to have no intrinsic length.

In heat transfer the half-line solution $w = 1 - \operatorname{erf}\bigl(x/2\sqrt{\kappa t}\bigr)$ is the textbook description of a semi-infinite body whose surface is suddenly brought to a new temperature; it gives the thermal penetration depth $\delta \simeq 2\sqrt{\kappa t}$ used in quenching, in welding and in estimating how deep the daily and annual temperature waves reach into soil and rock. The same formula, with $\kappa$ the chemical diffusivity, gives the carburising depth in steel and the junction depth in semiconductor dopant diffusion; because the profile depends only on $x/\sqrt{t}$ , doubling the depth costs four times the time, a scaling that is read directly off the similarity variable.

In geophysics and hydrology the Barenblatt solution describes the spreading of a finite mass of fluid through a porous medium — a groundwater mound, an oil or gas pocket, a spreading viscous gravity current — and its finite propagation speed, with a front at radius $\propto t^{\beta}$ , is what makes it usable: the linear diffusion equation would predict a trace of the pollutant everywhere immediately. Equations of porous-medium type also model radiative heat transport in a plasma, which is where Zeldovich and Kompaneetz met the solution.

In the theory of water waves the $t^{1/3}$ Airy similarity solution is the universal shape of the leading edge of a dispersive wave train: the oscillatory tail behind a front, with wavelength shrinking and amplitude decaying as $t^{-1/3}$ , is seen in ship wakes, in tsunami records, and in the collapse of an initial step for the Korteweg-de Vries equation (see travelling-wave solutions for the other half of that story, the solitons that run ahead of it).

Method

Travelling-wave solutions

Adapted

Add to PDF

Reduces a partial differential equation to an ordinary differential equation for a profile moving at constant speed, shows when the speed is free and when it is selected by the equation, derives the Korteweg-de Vries soliton and an explicit bistable reaction-diffusion front, and states the existence theorem for such fronts.

Prerequisites: First-order equations and characteristics, Flows on the line, Linearization and the Hartman-Grobman theorem

Definition(Travelling wave, plane wave, profile, speed)AdaptedAdapted from Evans (2010), §4.2.1, Eqs. (1)–(3), p. 167; Stein & Shakarchi (2003), Ch. 1, §1.1, p. 9

Let an equation be posed for $u(x,t)$ with $x \in \mathbb{R}$ . A solution of the form

u(x,t) = v(x - \sigma t) \qquad (x \in \mathbb{R},\ t \in \mathbb{R}) \tag{1}

is a travelling wave with speed $\sigma$ and profile $v$ . Its graph at any time is the graph of $v$ translated by $\sigma t$ : the shape never changes, only its position. With $x \in \mathbb{R}^n$ , a solution

u(x,t) = v(y\cdot x - \sigma t) , \qquad y \in \mathbb{R}^n\setminus\{0\} , \tag{2}

is a plane wave, with wavefronts the hyperplanes normal to $y$ and speed $\sigma/|y|$ in the direction $y/|y|$ . A travelling wave which decays at both ends, $v \to 0$ as $s \to \pm\infty$ , is a solitary wave; one which connects two distinct constant states, $v(-\infty) \ne v(+\infty)$ , is a front (or a kink). For a linear equation with constant coefficients it is convenient to allow complex-valued plane waves

u(x,t) = e^{i(y\cdot x - \sigma t)} ,

with $y$ the wave numbers and $\sigma$ the time frequency; substituting into the equation forces a relation $\sigma = \sigma(y)$ , the dispersion relation. Then $\sigma(y)/|y|$ is the phase velocity, the speed at which this single plane wave moves, and $D\sigma(y)$ is the group velocity, which is the speed at which a superposition $\int e^{i(y\cdot x-\sigma(y)t)}a(y)\,dy$ transports its energy: by stationary phase the main contribution at large $t$ with $x/t$ fixed comes from the wave numbers with $D\sigma(y) = x/t$ .

The simplest case is the linear wave equation $u_{tt} = u_{xx}$ : as Stein and Shakarchi observe, $F(x-t)$ and $F(x+t)$ solve it for every twice-differentiable $F$ , so both speeds $\pm1$ occur and every profile is admissible; and the change of variables $\xi = x+t$ , $\eta = x-t$ turns the equation into $v_{\xi\eta} = 0$ , whose double integration shows that every twice-differentiable solution is a sum $F(x+t) + G(x-t)$ of two such waves. The interesting questions arise when the equation restricts the profile, the speed, or both.

Method(The travelling-wave ansatz)AdaptedAdapted from Evans (2010), §4.2.1, pp. 169–171

Substitute. Put $u = v(s)$ with $s = x - \sigma t$ into the equation. Since $\partial_t u = -\sigma v'$ and $\partial_x^k u = v^{(k)}$ , every term becomes a term in $v$ and its derivatives: the result is an ordinary differential equation for $v$ , of the same order as the highest $x$ -derivative, containing $\sigma$ as a parameter.
Impose the conditions at infinity. A solitary wave asks $v, v', v'' \to 0$ as $s \to \pm\infty$ ; a front asks $v \to v_\pm$ as $s \to \pm\infty$ with $v' \to 0$ , where $v_\pm$ must be zeros of the reaction term, since they are constant solutions of the equation.
Reduce the order. Look for first integrals. A term $v^{(k)}$ integrates once immediately; a factor $v'$ makes the whole equation an exact derivative. For second-order autonomous problems, multiplying by $v'$ and integrating produces an energy relation of the form $\tfrac12(v')^2 + W(v) = \text{const}$ .
Use the conditions at infinity to kill the constants. This is where the ansatz earns its keep: each decay condition fixes one constant of integration, usually to zero, and what remains is a first-order separable equation.
Decide the status of $\sigma$ . Three outcomes occur, and they are not a matter of technique but of the equation:
- $\sigma$ is free: a travelling wave exists for every $\sigma$ in some range (the linear wave equation; the Korteweg-de Vries soliton);
- $\sigma$ is selected: exactly one value admits a solution with the required behaviour at $\pm\infty$ (the bistable front below);
- no travelling wave exists at all for the required boundary behaviour.
In the selected case the mechanism is a counting argument in the phase plane: the orbit leaving one equilibrium and the orbit entering the other each form a one-parameter family in $\sigma$ , and they coincide only for special $\sigma$ .
Verify, and read off the physics. Substitute back. Then extract what the answer says: the relation between amplitude and speed, the width of the front, and the direction of travel (the sign of $\sigma$ says which state advances).

Method(Existence of a front for a bistable nonlinearity)AdaptedAdapted from Evans (2010), §4.2.1c, hypotheses (13) and the conclusion at Eq. (26), pp. 170–175

Consider the scalar reaction-diffusion equation

u_t - u_{xx} = f(u) \qquad\text{in } \mathbb{R}\times(0,\infty) , \tag{3}

and suppose $f : \mathbb{R}\to\mathbb{R}$ is smooth and, for some $a$ with $0 < a < 1$ , satisfies

\text{(a) } f(0)=f(a)=f(1)=0; \quad \text{(b) } f<0 \text{ on } (0,a),\ f>0 \text{ on } (a,1); \quad \text{(c) } f'(0)<0,\ f'(1)<0; \quad \text{(d) } \int_0^1 f(z)\,dz>0 . \tag{4}

Then there is a speed $\sigma < 0$ and a profile $v$ for which $u(x,t) = v(x-\sigma t)$ solves (3) with

\lim_{s\to+\infty}v(s) = 1, \qquad \lim_{s\to-\infty}v(s) = 0, \qquad \lim_{s\to\pm\infty}v'(s) = 0 . \tag{5}

Evans adds that a more refined analysis, which he does not give, shows that this $\sigma$ is unique: for a nonlinearity satisfying (4) there is exactly one speed at which such a front exists. Contrast this with the soliton below, which exists for a whole half-line of speeds.

Idea of the proofOriginal

The ansatz turns (3) into the second-order equation $v'' + \sigma v' + f(v) = 0$ , which is the equation of a particle in the potential $-\int f$ with friction coefficient $\sigma$ . Written as a planar system it has saddle points at the two stable states $(0,0)$ and $(1,0)$ , because hypothesis (4)(c) makes the linearisation at each have eigenvalues of opposite sign. A front is exactly an orbit running from one saddle to the other, and the saddles each emit or absorb a single orbit, so the two curves must be made to coincide. The parameter $\sigma$ is the friction, and the energy $\tfrac12 w^2 + \int_0^v f$ changes along orbits at the rate $-\sigma w^2$ . With no friction the two curves miss each other in one direction, by hypothesis (4)(d); with enough friction they miss in the other; so somewhere in between they meet.

ProofAdaptedAdapted from Evans (2010), §4.2.1c, Eqs. (14)–(26), pp. 170–175

Evans presents this argument as an outline, without complete proofs of the geometric claims about the invariant curves; the steps below follow his route, and the places where a claim is asserted rather than established are flagged in the note after the proof.

Step 1 (the ansatz and the reduced equation). Putting $u = v(x-\sigma t)$ in (3) gives $-\sigma v' - v'' = f(v)$ , that is

v'' + \sigma v' + f(v) = 0 , \tag{6}

to be solved with the conditions (5). Note first why the limits must be $0$ and $1$ and not something else: a constant limit $v_\pm$ with $v'\to0$ forces $f(v_\pm) = 0$ , so the limits are among $0, a, 1$ ; and (4)(c) with $f'(a)\ge0$ makes $0$ and $1$ the stable constant solutions of (3) and $a$ the unstable one, so a front interpolating between stable states must join $0$ to $1$ .

Step 2 (phase plane). Set $w = v'$ . Then (6) becomes the autonomous planar system

v' = w, \qquad w' = -\sigma w - f(v) , \tag{7}

and the conditions (5) become $(v,w)\to(1,0)$ as $s\to+\infty$ and $(v,w)\to(0,0)$ as $s\to-\infty$ .

Step 3 (both equilibria are saddles). By (4)(a), $(0,0)$ and $(1,0)$ are the equilibria of (7) with $v \in \{0,1\}$ . Linearising (7) at $(v_*,0)$ gives the matrix $\begin{pmatrix}0 & 1\\ -f'(v_*) & -\sigma\end{pmatrix}$ , with characteristic equation $\lambda^2 + \sigma\lambda + f'(v_*) = 0$ and roots

\lambda^{\pm} = \frac{-\sigma \pm \bigl(\sigma^2 - 4f'(v_*)\bigr)^{1/2}}{2} .

By hypothesis (4)(c), $f'(v_*)<0$ at both $v_*=0$ and $v_*=1$ , so the discriminant exceeds $\sigma^2 \ge 0$ and the roots are real, while their product $f'(v_*)$ is negative, so they have opposite signs. Hence both equilibria are saddle points (see linearisation). Consequently exactly one orbit leaves $(0,0)$ into the quadrant $\{v>0,w>0\}$ — call its closure the unstable curve $W^{\mathrm{u}}$ , tangent at $(0,0)$ to the line $w = \lambda^+_0 v$ — and exactly one orbit enters $(1,0)$ from that quadrant — the stable curve $W^{\mathrm{s}}$ , tangent at $(1,0)$ to $w = \lambda^-_1(v-1)$ . Both curves and both eigenvalues depend on $\sigma$ .

Step 4 (what is to be shown). A front exists precisely when, for some $\sigma$ ,

W^{\mathrm{u}} = W^{\mathrm{s}} \quad\text{in } \{v>0,\ w>0\} ,

for then that common orbit is a heteroclinic orbit from $(0,0)$ to $(1,0)$ , and reading it as a graph $w = w(v)$ and integrating $v' = w$ recovers $v$ .

Step 5 (an energy that friction drives). Define

E(v,w) = \frac{w^2}{2} + \int_0^v f(z)\,dz .

Along a solution of (7),

\frac{d}{ds}E\bigl(v(s),w(s)\bigr) = w\,w' + f(v)\,v' = w\bigl(-\sigma w - f(v)\bigr) + f(v)\,w = -\sigma\,w^2 .

So for $\sigma<0$ , $E$ is non-decreasing along every orbit, and for $\sigma = 0$ it is constant, so the orbits lie in the level sets of $E$ .

Step 6 (both curves cross a fixed vertical line). Fix a small $\varepsilon>0$ with $a+\varepsilon<1$ and let $L$ be the vertical segment above $(a+\varepsilon,0)$ . Let $R$ be the region whose three straight sides are the segment of $w=0$ , the segment of $v=0$ , and $L$ . Evans argues, for $\sigma<0$ , that $W^{\mathrm{u}}$ enters $R$ from the origin and can leave it only across $L$ , at some height $w_0(\sigma)$ , and that $W^{\mathrm{s}}$ likewise reaches $L$ , at some height $w_1(\sigma)$ .

Step 7 (at zero friction the stable curve is higher). For $\sigma = 0$ each curve lies in a level set of $E$ by Step 5. On $W^{\mathrm{u}}$ that level is $E(0,0) = 0$ , so at $v = a+\varepsilon$ ,

\frac{w_0(0)^2}{2} = -\int_0^{a+\varepsilon} f , \qquad\text{hence}\qquad w_0(0) = \Bigl(-2\int_0^{a+\varepsilon} f\Bigr)^{1/2} ,

which is a positive real number because $f<0$ on $(0,a)$ and $\varepsilon$ is small. On $W^{\mathrm{s}}$ the level is $E(1,0) = \int_0^1 f$ , so

w_1(0) = \Bigl(2\int_0^1 f - 2\int_0^{a+\varepsilon} f\Bigr)^{1/2} .

Since $\int_0^1 f>0$ by hypothesis (4)(d), the quantity under the second root exceeds that under the first, and therefore

w_0(0) < w_1(0) . \tag{8}

This is the only place (4)(d) is used, and it is what gives the front its direction.

Step 8 (at large friction the unstable curve is higher). Fix $\beta>0$ and consider the segment $T = \{0\le v\le a+\varepsilon,\ w = \beta v\}$ . On $T$ , dividing the two equations of (7),

\frac{w'}{v'} = \frac{-\sigma w - f(v)}{w} = -\sigma - \frac{f(v)}{\beta v} .

Because $f$ is smooth with $f(0)=0$ , the ratio $|f(v)/v|$ is bounded by some $C$ on $0\le v\le a+\varepsilon$ , so $w'/v' \ge -\sigma - C/\beta$ , which exceeds $\beta$ once $\sigma<0$ with $|\sigma|$ large enough. The orbit therefore crosses $T$ more steeply than $T$ itself and cannot escape below it, so $W^{\mathrm{u}}$ stays above the line $w = \beta v$ until it reaches $L$ , giving $w_0(\sigma) \ge \beta(a+\varepsilon)$ . On the other hand Evans asserts $w_1(\sigma) \le w_1(0)$ for all $\sigma\le0$ . Choosing $\beta$ large enough that $\beta(a+\varepsilon) > w_1(0)$ and then $\sigma$ sufficiently negative,

w_0(\sigma) > w_1(\sigma) . \tag{9}

Step 9 (intermediate value). The two heights $w_0$ and $w_1$ depend smoothly on $\sigma$ , so by (8), (9) and continuity there is a $\sigma<0$ with

w_0(\sigma) = w_1(\sigma) .

At that $\sigma$ the two curves meet at a point of $L$ ; since both are orbits of the same autonomous system (7), and an orbit through a point is unique, they are the same orbit. This is the heteroclinic connection of Step 4, and it yields a travelling wave of (3) satisfying (5). $\square$

NoteOriginal

Three claims in the proof are asserted by Evans rather than proved, and are repeated here on his authority: that $W^{\mathrm{u}}$ cannot leave the region $R$ except through $L$ and that $W^{\mathrm{s}}$ likewise reaches $L$ (Step 6); that $w_1(\sigma)\le w_1(0)$ for every $\sigma\le0$ (Step 8); and that $w_0$ and $w_1$ depend smoothly on $\sigma$ (Step 9). He also states the uniqueness of the selected speed as the outcome of a more refined analysis which he does not present. The explicit example below is independent of all of this: for one particular cubic $f$ satisfying (4), it exhibits the front and its speed in closed form and verifies them by substitution.

Remark(Plane waves and dispersion)AdaptedAdapted from Evans (2010), §4.2.1a, pp. 167–169

For linear constant-coefficient equations, inserting $u = e^{i(y\cdot x - \sigma t)}$ is step 1 of the method with an exponential profile, and the resulting dispersion relation already tells one most of what the equation does.

Equation	Dispersion relation	What it says
$u_t - \Delta u = 0$	$\sigma = -i\lvert y\rvert^2$	purely imaginary $\sigma$ : the factor $e^{-\lvert y\rvert^2t}$ damps, and damps faster the shorter the wave
$u_{tt} - \Delta u = 0$	$\sigma = \pm\lvert y\rvert$	real $\sigma$ , with speed $\lvert\sigma\rvert/\lvert y\rvert = 1$ for every $y$ : no damping, no dispersion
$u_{tt} - \Delta u + m^2u = 0$	$\sigma = \pm(\lvert y\rvert^2+m^2)^{1/2}$	real, but the speed now depends on $\lvert y\rvert$ : dispersion, with the slowly oscillating waves travelling faster
$iu_t + \Delta u = 0$	$\sigma = \lvert y\rvert^2$	dispersion
$u_t + u_{xxx} = 0$	$\sigma = -y^3$	dispersion, with speed $y^2$

Damping shows up as $\operatorname{Im}\sigma \ne 0$ and dispersion as a non-constant phase velocity $\sigma(y)/\lvert y\rvert$ . Dispersion is what makes the nonlinear balance of the next example possible: in the Korteweg-de Vries equation the dispersive term $u_{xxx}$ spreads a pulse out while the nonlinear term $uu_x$ steepens it, and a soliton is the shape for which the two exactly cancel.

Example(The Korteweg-de Vries soliton)AdaptedAdapted from Evans (2010), §4.2.1b, Eqs. (4)–(11), pp. 169–170; Strauss (2008), §14.2, Eqs. (1)–(3), pp. 390–391

Take the Korteweg-de Vries equation in Evans's and Strauss's normalisation,

u_t + 6uu_x + u_{xxx} = 0 \qquad\text{in } \mathbb{R}\times(0,\infty) ,

a model for long surface waves on shallow water, and look for a solitary wave $u = v(x-\sigma t)$ .

Substitute. With $s = x - \sigma t$ , the equation becomes $-\sigma v' + 6vv' + v''' = 0$ .
Integrate once. Every term is a derivative: $\bigl(-\sigma v + 3v^2 + v''\bigr)' = 0$ , so
$-\sigma v + 3v^2 + v'' = a$
for a constant $a$ .
Integrate again. Multiply by $v'$ . Then $-\sigma vv' + 3v^2v' + v''v' = av'$ , and each side is a derivative, giving
$\frac{(v')^2}{2} = -v^3 + \frac{\sigma}{2}v^2 + av + b$
for a second constant $b$ .
Impose decay. A solitary wave has $v, v', v'' \to 0$ as $s \to \pm\infty$ . Letting $s\to\infty$ in step 2 gives $a = 0$ , and then in step 3 gives $b = 0$ . Hence
$(v')^2 = v^2\bigl(\sigma - 2v\bigr) ,$
which already forces $\sigma > 0$ for a non-trivial non-negative profile, and bounds the amplitude by $\sigma/2$ .
Solve. Separating variables and substituting $v = \tfrac{\sigma}{2}\operatorname{sech}^2\theta$ reduces the integral to $ds = (2/\sqrt{\sigma})\,d\theta$ , so $\theta = \tfrac{\sqrt\sigma}{2}(s-c)$ and
$u(x,t) = \frac{\sigma}{2}\,\operatorname{sech}^2\!\left(\frac{\sqrt\sigma}{2}\,\bigl(x - \sigma t - c\bigr)\right) \qquad (\sigma>0,\ c \in \mathbb{R}) .$
Substituting this expression into the equation confirms it identically.
What the formula says. There is a soliton for every speed $\sigma>0$ : the speed is free, the opposite of the bistable case. Its amplitude is $\sigma/2$ and its width is of order $2/\sqrt\sigma$ , so taller solitons are narrower and faster, as Strauss puts it: tall, thin and fast for large $\sigma$ ; short, fat and slow for small $\sigma$ . Numerically, $\sigma = 2$ gives amplitude $1$ and $\sigma = 8$ gives amplitude $4$ at twice the width-inverse, so a soliton four times as tall travels four times as fast and is half as wide.
Conventions differ, and the formula with them. Arfken, Weber and Harris write the equation as $\psi_t + \psi\psi_x + \psi_{xxx} = 0$ , without the factor $6$ , and obtain $\psi = 3c\,\operatorname{sech}^2\bigl(\tfrac12\sqrt{c}\,(x-ct)\bigr)$ — amplitude $3c$ rather than $c/2$ . The two are the same family: replacing $u$ by $\psi/6$ turns one equation into the other, and $\tfrac{\sigma}{2}\mapsto 3\sigma$ under that substitution. Strauss remarks that the factor $6$ has no special significance; it is chosen to tidy the formula. Always check which normalisation a quoted soliton formula belongs to.
Why solitons are remarkable. Strauss records the history: J. S. Russell saw a solitary heap of water run along a canal for a mile or two without changing shape in 1834, and Kruskal and Zabusky found numerically in the 1960s that when a fast soliton overtakes a slow one the two emerge from the nonlinear interaction with their shapes and speeds intact, delayed only by a phase shift, and that general initial data decompose as $t\to\infty$ into finitely many solitons plus a decaying dispersive tail. The explanation is that the equation is integrable, with an exact solution available for general data through inverse scattering; the three invariants $\int u$ , $\int u^2$ and $\int(\tfrac12u_x^2 - u^3)$ are the first of infinitely many conserved quantities.

ExampleOriginal

An explicit bistable front, with its selected speed. Take the cubic

f(u) = u(1-u)(u-a), \qquad 0 < a < \tfrac12 ,

in equation (3). It satisfies all four hypotheses (4): the zeros are $0, a, 1$ ; $f<0$ on $(0,a)$ and $f>0$ on $(a,1)$ by inspection of the three factors; $f'(0) = -a < 0$ and $f'(1) = -(1-a) < 0$ ; and

\int_0^1 f(z)\,dz = \int_0^1 \bigl(-z^3 + (1+a)z^2 - az\bigr)\,dz = -\frac14 + \frac{1+a}{3} - \frac{a}{2} = \frac{1}{12} - \frac{a}{6} ,

which is positive exactly when $a < \tfrac12$ . So the theorem applies and predicts a unique $\sigma<0$ .

Guess the shape. The reduced equation is $v'' + \sigma v' + v(1-v)(v-a) = 0$ with $v(-\infty)=0$ , $v(+\infty)=1$ . Try an orbit on which $w = v'$ is the quadratic $v' = \kappa\,v(1-v)$ , with $\kappa>0$ : this vanishes at both equilibria and is positive between them, so it is a candidate for the heteroclinic orbit of Step 4.
Substitute. Differentiating, $v'' = \kappa(1-2v)v' = \kappa^2 v(1-v)(1-2v)$ . Putting both into the reduced equation,
$\kappa^2 v(1-v)(1-2v) + \sigma\kappa\,v(1-v) + v(1-v)(v-a) = 0 ,$
and dividing by $v(1-v)$ , which is nonzero for $0<v<1$ ,
$\kappa^2(1-2v) + \sigma\kappa + (v-a) = 0 .$
Match powers of $v$ . The coefficient of $v$ is $1 - 2\kappa^2$ and the constant term is $\kappa^2 + \sigma\kappa - a$ . Both must vanish, so
$\kappa = \frac{1}{\sqrt2}, \qquad \sigma = \frac{a - \kappa^2}{\kappa} = \sqrt2\Bigl(a - \tfrac12\Bigr) = \frac{2a-1}{\sqrt2} .$
The guess has worked, and it has produced one value of $\sigma$ and no freedom: the speed is selected, exactly as the theorem says. And $\sigma<0$ precisely when $a<\tfrac12$ , which is precisely when $\int_0^1f>0$ — the hypothesis used in Step 7 of the proof.
Integrate for the profile. Solving $v' = v(1-v)/\sqrt2$ by separation gives the logistic profile
$v(s) = \frac{1}{1 + e^{-s/\sqrt2}} , \qquad\text{so}\qquad u(x,t) = \frac{1}{1 + \exp\bigl(-(x-\sigma t)/\sqrt2\bigr)} ,$
with $v(-\infty)=0$ , $v(+\infty)=1$ , $v(0)=\tfrac12$ and $v'(\pm\infty)=0$ , as (5) requires. Equivalently $v(s) = \tfrac12\bigl[1+\tanh\bigl(s/(2\sqrt2)\bigr)\bigr]$ .
Check. Substituting this $u$ into $u_t - u_{xx} = u(1-u)(u-a)$ with $\sigma = (2a-1)/\sqrt2$ gives zero identically, for every $a$ .
Numbers. For $a = \tfrac14$ : $\sigma = -\tfrac{1}{2\sqrt2} = -0.3535534$ , so the front moves to the left at that speed, meaning the state $u=1$ advances into the state $u=0$ — which is the state with the larger value of $\int f$ , as it should be. The front is centred where $u = \tfrac12$ , at $x = \sigma t$ , and its width is of order $\sqrt2$ : $u = 0.8044$ at $s = 2$ and $u = 0.1956$ at $s = -2$ . At the balanced value $a=\tfrac12$ the speed is zero and the front stands still; for $a>\tfrac12$ it reverses, and $u=0$ invades instead.
A caution on generality. This is one explicit member of the family the theorem covers, not a proof of the theorem: the trick of guessing $v' = \kappa v(1-v)$ works because the nonlinearity is exactly cubic with zeros at $0$ and $1$ , and it has no analogue for general $f$ . Note also that if $f$ is positive on all of $(0,1)$ , so that $f'(0)>0$ and hypothesis (4)(c) fails, none of the above applies: the equilibrium at the origin is then not a saddle, Step 3 breaks down, and the speed is no longer determined by this argument.

ApplicationOriginal

Fronts of exactly the bistable kind are how a transition spreads through a medium that has two stable states. In combustion, $u$ is a scaled temperature or reactant fraction and the front is a premixed flame, whose selected speed is the laminar burning velocity — a material property, measured in cm/s, that the equation predicts from the reaction kinetics and the thermal diffusivity. In materials science the same equation with $f = -W'(u)$ for a double-well $W$ is the Allen-Cahn equation, and the front is the interface between two phases of a solidifying alloy or two orientations of a magnetic domain; the unequal well depths, which are hypothesis (4)(d), drive the interface towards the deeper well at the speed computed above, and the front width $\sqrt{2}$ in the example is the interface thickness in units of the diffusion length. In neurophysiology the cubic $f$ is the fast part of the FitzHugh-Nagumo reduction of the Hodgkin-Huxley equations, and the front is the leading edge of an action potential travelling along an axon; that the conduction velocity is a definite number, not a free parameter, is exactly the speed selection of the statement above. In population biology the same structure describes an advantageous gene or an invading species spreading into territory held by a resident one.

Solitary waves appear wherever weak nonlinearity balances weak dispersion. The Korteweg-de Vries soliton is the shallow-water wave Russell chased on horseback, and it is also the model for internal waves in a stratified ocean and for ion-acoustic pulses in a plasma; the closely analogous soliton of the nonlinear Schrödinger equation is the basis of soliton transmission in optical fibres, where the pulse shape survives propagation because self-phase modulation cancels group-velocity dispersion. Dispersion relations themselves, the linear part of this unit, are the daily tool for reading a wave-bearing equation: whether it damps, at what rate the short waves go, and at what speed a wave packet actually carries energy, which is the group velocity and not the phase velocity.