Method of Characteristics for First-Order PDEs

The Method of Characteristics converts a first-order partial differential equation (PDE) into a coupled system of ordinary differential equations (ODEs) along special curves called characteristics, along which the PDE collapses to a simple total-derivative relation. Key result: the solution surface of $a u_{x} + b u_{y} = c$ is foliated by curves on which information propagates, reducing the problem to integrating ODEs. Applies to quasilinear and fully nonlinear first-order PDEs in any number of variables. Origin: developed by Cauchy, Monge, and Charpit in the late 18th–19th centuries.

The Method of Characteristics is a technique for solving first-order PDEs by identifying curves in the domain—characteristic curves—along which the PDE becomes an ordinary differential equation. For a quasilinear equation $a (x, y, u) u_{x} + b (x, y, u) u_{y} = c (x, y, u)$ , the directional derivative of $u$ in the direction $(a, b)$ equals $c$ ; this geometric observation transforms the PDE into a system of ODEs that can be integrated step by step. The method matters because first-order PDEs model transport, wave propagation, traffic flow, and optics, and characteristics reveal how information travels through these systems. A surprising fact: even for perfectly smooth initial data, characteristics of nonlinear equations can cross, forcing solutions to develop discontinuities (shocks) in finite time—the breakdown of classical solutions is built into the geometry.

The Core Insight

A first-order PDE prescribes the rate of change of $u$ in one special direction at each point. If you move along a curve tangent to that direction, the PDE tells you exactly how $u$ changes—an ODE. By stitching together these one-dimensional solutions, you reconstruct the entire solution surface. The PDE is really a family of ODEs in disguise.

Formal Statement

Quasilinear First-Order PDE

A quasilinear first-order PDE in two independent variables has the form

a (x, y, u) u_{x} + b (x, y, u) u_{y} = c (x, y, u)

⟳This is the most general first-order PDE that is linear in the highest (first) derivatives

u_{x}, u_{y}

; coefficients may depend nonlinearly on

u

itself.

Here $u = u (x, y)$ is the unknown, $u_{x} = \partial u / \partial x$ , $u_{y} = \partial u / \partial y$ , and $a, b, c$ are given functions. The solution is a surface $z = u (x, y)$ in $(x, y, z)$ -space.

Characteristic Curves and the Characteristic System

The characteristic curves are the projections onto the $(x, y)$ -plane of the integral curves of the characteristic system of ODEs:

\frac{d x}{d t} = a (x, y, u), \frac{d y}{d t} = b (x, y, u), \frac{d u}{d t} = c (x, y, u)

⟳Parametrize a curve

(x (t), y (t))

in the domain; the chain rule gives

d u / d t = u_{x} \overset{x}{˙} + u_{y} \overset{y}{˙}

. Choosing

\overset{x}{˙} = a

\overset{y}{˙} = b

makes the right-hand side equal

a u_{x} + b u_{y} = c

The full three-dimensional curves $(x (t), y (t), u (t))$ lying on the solution surface are called characteristic strips (in the fully nonlinear case) or simply characteristics.

Cauchy Problem and Initial Curve

The Cauchy problem supplies initial data on a curve $Γ$ in the domain, parametrized by $s$ :

x = x_{0} (s), y = y_{0} (s), u = u_{0} (s) at t = 0

⟳We prescribe the value of

u

along a one-dimensional curve

Γ

; characteristics emanate from each point of

Γ

, sweeping out a two-dimensional surface as

(s, t)

vary.

The solution surface is generated by launching one characteristic from each point of $Γ$ .

Existence and Uniqueness via the Transversality Condition

Let $a, b, c \in C^{1}$ and let the initial curve $Γ$ be $C^{1}$ . If the transversality (non-characteristic) condition

J = det (x_{0}^{'} (s) y_{0}^{'} (s) a (x_{0}, y_{0}, u_{0}) b (x_{0}, y_{0}, u_{0})) = x_{0}^{'} (s) b - y_{0}^{'} (s) a \neq = 0

⟳

J

is the Jacobian

\partial (x, y) / \partial (s, t)

t = 0

; it is nonzero precisely when

Γ

is nowhere tangent to a characteristic direction, so the map

(s, t) \mapsto (x, y)

is locally invertible by the inverse function theorem.

holds along $Γ$ , then there exists a unique $C^{1}$ solution $u (x, y)$ in a neighborhood of $Γ$ satisfying the PDE and the initial data.

If $J = 0$ on $Γ$ —i.e., the initial curve is a characteristic—the Cauchy problem has either no solution or infinitely many. Data must be consistent with the ODE $d u / d t = c$ for any solution to exist at all. Always check the transversality condition first.

Why It Holds / Derivation

Reduction of the PDE to ODEs

Step 1 — Geometry of the solution surface. Suppose

z = u (x, y)

solves [1]. Write the surface implicitly as

F (x, y, z) = u (x, y) - z = 0

. Its normal vector is

\nabla F = (u_{x}, u_{y}, - 1) .

⟳Take gradients of

F = u (x, y) - z

; the

z

-component is

\partial (- z) / \partial z = - 1

Step 2 — The PDE as orthogonality. Equation [1] states

a u_{x} + b u_{y} - c = 0

, which is exactly

(a, b, c) \cdot (u_{x}, u_{y}, - 1) = 0.

⟳Expand the dot product; the third term is

c \cdot (- 1) = - c

, and setting the sum to zero recovers [1].

So the vector field $(a, b, c)$ is tangent to the solution surface at every point.

Step 3 — Integral curves lie on the surface. A curve whose tangent is always

(a, b, c)

—i.e., a solution of system [2]—stays within the surface once it starts there, because its velocity is perpendicular to the surface normal:

\frac{d}{d t} F (x (t), y (t), u (t)) = \nabla F \cdot (\overset{x}{˙}, \overset{y}{˙}, \overset{u}{˙}) = (u_{x}, u_{y}, - 1) \cdot (a, b, c) = 0.

⟳Apply the chain rule to

F

along the curve, then substitute [2] for the velocity and use [6].

Hence $F$ is constant ( $= 0$ ) along characteristics: they remain on the solution surface.

Step 4 — Building the surface. Launch a characteristic from each point

(x_{0} (s), y_{0} (s), u_{0} (s))

Γ

. This yields

x = x (s, t)

y = y (s, t)

u = u (s, t)

. When

J \neq = 0

(condition [4]), invert to get

s = s (x, y)

t = t (x, y)

and write

u (x, y) = u (s (x, y), t (x, y))

. This function is

C^{1}

and satisfies [1].

■

Augustin-Louis Cauchy | 1789–1857

Cauchy formalized the existence theory for first-order PDEs and laid the groundwork for what became the Cauchy problem. Building on Gaspard Monge's geometric picture of PDEs as families of surfaces and Paul Charpit's treatment of the fully nonlinear case (the Lagrange–Charpit method), Cauchy's rigorous approach—using the propagation of data along characteristic curves—turned a collection of geometric tricks into a systematic theory that underlies modern hyperbolic PDE analysis.

Building Intuition

Imagine dye dropped into a flowing river. The fluid velocity field is $(a, b)$ , and the dye concentration $u$ rides along with the current. A characteristic curve is simply a path traced by a floating particle. The equation $d u / d t = c$ tells you how the concentration changes as you drift downstream—if $c = 0$ , the dye value is carried unchanged along the flow (pure transport). The solution everywhere is just "look back upstream along your characteristic to where you started, and read off the initial value (adjusted by integrating $c$ )." Information cannot cross characteristics; each one is an isolated channel of communication.

Characteristics emanating from an initial curve $\Gamma$; the solution at point $P$ is determined by tracing back along its characteristic to $\Gamma$.

Worked Examples

Linear Transport Equation

Given: Solve

u_{t} + 2 u_{x} = 0

with initial condition

u (x, 0) = sin x

. (Here

y = t

a = 1

in the

t

-direction,

b = 2

c = 0

.) Step 1 — Characteristic ODEs. Treating

t

as the parameter:

\frac{d t}{d τ} = 1, \frac{d x}{d τ} = 2, \frac{d u}{d τ} = 0.

⟳Match

u_{t} + 2 u_{x} = 0

a u_{t} + b u_{x} = c

with

a = 1, b = 2, c = 0

, then read off [2].

Step 2 — Integrate. From

d t / d τ = 1

take

t = τ

. Then

x = 2 t + x_{0}

and

u = u_{0}

(constant along characteristics). Step 3 — Apply data. At

t = 0

x = x_{0}

and

u_{0} = sin x_{0}

. Invert:

x_{0} = x - 2 t

. Step 4 — Solution.

u (x, t) = sin (x - 2 t) .

⟳Substitute

x_{0} = x - 2 t

into

u = u_{0} = sin x_{0}

; the profile

sin

translates rightward at speed 2.

u (x, t) = sin (x - 2 t) (a wave moving right at speed 2 units/time)

Quasilinear Equation with a Source

Given: Solve

u_{x} + u_{y} = u

with

u (x, 0) = x

along

Γ : (x_{0} (s), y_{0} (s), u_{0} (s)) = (s, 0, s)

. Step 1 — Check transversality.

J = x_{0}^{'} b - y_{0}^{'} a = (1) (1) - (0) (1) = 1 \neq = 0

. ✓ Non-characteristic. Step 2 — Characteristic system.

\frac{d x}{d t} = 1, \frac{d y}{d t} = 1, \frac{d u}{d t} = u .

⟳Read

a = 1, b = 1, c = u

directly from

u_{x} + u_{y} = u

and insert into [2].

Step 3 — Integrate with initial data.

x = s + t, y = t, u = s e^{t} .

⟳

x (0) = s, y (0) = 0, u (0) = s

; the equation

d u / d t = u

gives exponential growth

u = u (0) e^{t}

Step 4 — Invert.

t = y

and

s = x - t = x - y

. Step 5 — Solution.

u (x, y) = (x - y) e^{y} .

⟳Substitute

s = x - y

t = y

into

u = s e^{t}

u (x, y) = (x - y) e^{y}

Check:

u_{x} = e^{y}

u_{y} = - e^{y} + (x - y) e^{y}

; sum

= (x - y) e^{y} = u

. ✓

Common Misconceptions

"Characteristics are always straight lines." False. They are straight only when

a

and

b

are constant (or, in the quasilinear case, when

u

stays constant along them so the slope

b / a

does not change). For genuinely nonlinear equations like the inviscid Burgers' equation

u_{t} + u u_{x} = 0

, characteristic slopes depend on

u

, so they generically bend, converge, and intersect.

Shock Formation in Burgers' Equation

Consider $u_{t} + u u_{x} = 0$ with $u (x, 0) = - x$ on $[- 1, 1]$ . Characteristics satisfy $d x / d t = u = u_{0} = - x_{0}$ , so

x = x_{0} - x_{0} t = x_{0} (1 - t), u = - x_{0} .

⟳

u

is constant along characteristics, so its initial value

- x_{0}

persists; the speed equals that constant.

At $t = 1$ all characteristics collapse to $x = 0$ : the map $(x_{0}, t) \mapsto x$ becomes non-invertible ( $J = 1 - t = 0$ ). The classical solution ceases to exist at $t = 1$ even though the initial data was smooth. Beyond this, a weak solution with a shock wave is required, governed by the Rankine–Hugoniot condition. This shows the method's local nature: smooth solutions need not extend globally.

Connections

The Method of Characteristics is the gateway to several neighboring areas. It generalizes the Linear Transport Equation to nonlinear settings and underpins the theory of Conservation Laws, where shock formation and weak solutions take center stage. The fully nonlinear version connects to Hamilton-Jacobi Equations via the Lagrange–Charpit equations, which mirror Hamilton's canonical equations in Classical Mechanics. Geometrically it ties to First-Order ODE Systems (the characteristic system itself) and to the Inverse Function Theorem (the transversality condition). In numerical analysis, it informs upwind schemes and the CFL Condition for stability.

Real-World Applications

Traffic Flow (Lighthill–Whitham–Richards model). The conservation law $ρ_{t} + (ρ v (ρ))_{x} = 0$ for car density $ρ$ is solved via characteristics; shock formation models the onset of traffic jams.
Gas Dynamics. Compressible flow and supersonic nozzle design use characteristic meshes (method of characteristics for the 2D steady Euler equations) to compute shock-free contours.
Geometrical Optics. Light rays are the characteristics of the eikonal equation $∣\nabla S ∣^{2} = n^{2}$ , a Hamilton–Jacobi PDE; this connects wave fronts to ray tracing.
Population Dynamics (McKendrick–von Foerster equation). Age-structured population PDEs $n_{t} + n_{a} = - μ n$ are integrated exactly along characteristics of slope 1 in the age–time plane.
Chromatography & Reactive Transport. Solute transport with adsorption is modeled by quasilinear PDEs whose characteristics predict band velocities and elution profiles.
Financial Mathematics. Certain first-order Hamilton–Jacobi–Bellman equations in optimal control reduce to characteristic ODEs along optimal trajectories.
Computer Graphics. The eikonal/level-set framework for distance fields and shape propagation relies on characteristics (fast marching methods).

Historical Development

1770s — Joseph-Louis Lagrange develops methods for linear and quasilinear first-order PDEs, introducing the idea of auxiliary ODEs (Lagrange's method).
1784 — Paul Charpit de Villecourt submits a memoir treating the general (fully nonlinear) first-order PDE; though never published, it founds the Lagrange–Charpit method.
1770s–1810s — Gaspard Monge gives the geometric interpretation of PDEs via families of surfaces, characteristic curves, and the "Monge cone."
1819 — Pierre-Simon Laplace and contemporaries refine the connection to envelope and contact geometry.
1830s — Augustin-Louis Cauchy establishes rigorous existence theory for the Cauchy problem along non-characteristic data.
1860s — William Rowan Hamilton's earlier optics work (1820s–30s) is recognized as the characteristic theory of the eikonal/Hamilton–Jacobi equation.
1900s — Sophus Lie embeds the theory within contact transformations and Lie group symmetry.
1950s — Peter Lax, Olga Oleinik, and others develop weak solutions, entropy conditions, and shock theory for nonlinear conservation laws, extending characteristics past their breakdown.
1980s — Crandall & Lions introduce viscosity solutions, giving a rigorous global theory for Hamilton–Jacobi equations where classical characteristics fail.

Subtleties & Limits

The Method of Characteristics is fundamentally local: it guarantees a smooth solution only in a neighborhood of the initial curve. For nonlinear equations, characteristics may cross in finite time, producing shocks, and the classical solution no longer exists. The modern resolution—weak solutions with entropy/viscosity conditions—restores global existence and uniqueness, but identifying the physically correct weak solution (entropy condition selection) and extending the theory to systems of conservation laws in multiple space dimensions remains deep and partially open. Global well-posedness for multidimensional systems of hyperbolic conservation laws (e.g., compressible Euler with large data) is a major unsolved problem. Likewise, characteristics for fully nonlinear higher-dimensional PDEs intertwine with the still-active theory of viscosity solutions and weak KAM theory.