Ordinary differential equations

Lens

MotivationH

Laws of change

Most laws of nature do not tell you what a quantity is; they tell you how fast it changes. Radioactive carbon decays at a rate proportional to the amount present — which is how archaeologists date charcoal from Mungo in Australia or bone from Olduvai Gorge. A cup of chai cools at a rate proportional to its excess temperature over the room. An epidemic spreads at a rate proportional to the number of contacts between infected and susceptible people. A pendulum's acceleration is proportional to the sine of its angle.

Each law is an ordinary differential equation (ODE): an equation relating an unknown function of one variable (usually time) to its derivatives. Solving it means finding the whole function. Three kinds of question arise:

Explicit solution. Some equations can be solved by formulas; this lesson covers the main classes.
Existence and uniqueness. Most equations cannot be solved by formulas. Is there a solution anyway? Only one? For how long?
Qualitative behaviour. Without a formula, what can we say — does the solution settle down, oscillate, blow up?

This lesson assumes calculus of one variable and linear algebra (eigenvalues and eigenvectors).

IntuitionH

Following the arrows

An equation $x' = f(t, x)$ assigns a slope $f(t, x)$ to every point of the $(t, x)$ -plane. Draw a short segment with that slope at many points: this is a slope field. A solution is a curve that is tangent to the segments everywhere — you "follow the arrows". Euler's method makes this literal: take a small step $h$ in the direction of the arrow, look again, repeat.

The logistic equation x′ = x(1 − x)

Click to draw solutions. Every solution starting with 0 < x < 1 rises towards 1; none ever crosses the lines x = 0 or x = 1.

Examples

dy/dx = y·(1 − y)

Click the field to draw a solution, or type a start point below.

Start pointx₀y₀

Solution curves: 4. They are numerical approximations (fourth-order Runge–Kutta), not exact formulas.

Latest curve through (0, −0.05).

Notice two things. Solutions never cross: if two solutions met at a point, they would be two solutions of the same initial value problem, and the uniqueness theorem below forbids that (for well-behaved $f$ ). And the horizontal lines $x = 0$ and $x = 1$ are solutions themselves — equilibria — which the other solutions approach or flee. Much of the theory is visible in the picture before any formula is written.

DefinitionsH

Definitions

Definition— ordinary differential equation

An ordinary differential equation of order $n$ is an equation $F\big(t, x, x', \dots, x^{(n)}\big) = 0$ relating an unknown function $x$ of one variable $t$ to its first $n$ derivatives. A solution on an interval $I$ is an $n$ times differentiable function $x : I \to \R$ (or $\R^d$ ) satisfying the equation for every $t \in I$ . An initial value problem (IVP) adds conditions $x(t_0) = x_0, \dots, x^{(n-1)}(t_0) = x_{n-1}$ .

The equation is linear if it has the form $a_n(t)x^{(n)} + \cdots + a_1(t)x' + a_0(t)x = g(t)$ , and homogeneous if moreover $g = 0$ . A first-order system is $x' = f(t, x)$ with $x(t) \in \R^d$ ; it is autonomous if $f$ does not depend on $t$ . Every $n$ -th order equation $x^{(n)} = G(t, x, \dots, x^{(n-1)})$ becomes a first-order system for the vector $(x, x', \dots, x^{(n-1)})$ — for example, $x'' = -x$ becomes $(x, v)' = (v, -x)$ — so first-order systems are the general case.

Definition— equilibrium and stability

An equilibrium of the autonomous system $x' = f(x)$ is a point $x^*$ with $f(x^*) = 0$ ; the constant function $x(t) = x^*$ is then a solution. It is stable if solutions starting close to $x^*$ stay close for all $t \ge 0$ (for every $\varepsilon > 0$ there is $\delta > 0$ such that $|x(0) - x^*| < \delta$ implies $|x(t) - x^*| < \varepsilon$ for all $t \ge 0$ ), and asymptotically stable if moreover such solutions tend to $x^*$ .

Definition— Lipschitz condition

$f(t, x)$ is Lipschitz in $x$ on a set $R$ , with constant $K$ , if $|f(t, x) - f(t, y)| \le K|x - y|$ for all $(t, x), (t, y) \in R$ . By the mean value theorem, this holds on any closed box where $\partial f/\partial x$ is continuous (take $K$ = the maximum of its absolute value, or of the norm of the Jacobian matrix in several variables).

ProofH

Equations we can solve

Separable equations $x' = g(t)\,h(x)$ . Wherever $h(x) \ne 0$ we may write $\frac{x'}{h(x)} = g(t)$ . If $H' = 1/h$ and $G' = g$ , the chain rule shows that $\frac{d}{dt}H(x(t)) = \frac{x'(t)}{h(x(t))} = g(t)$ , so $H(x(t)) = G(t) + C$ ; solve for $x$ . In addition every zero $x^*$ of $h$ gives a constant solution $x \equiv x^*$ . (The familiar shorthand " $\int dx/h(x) = \int g(t)\,dt$ " is exactly this argument.)

Example— Three separable equations, three behaviours

Logistic growth $x' = x(1 - x)$ . Partial fractions $\frac{1}{x(1-x)} = \frac1x + \frac1{1-x}$ give $\ln\left|\frac{x}{1 - x}\right| = t + c$ , so $x(t) = \dfrac{1}{1 + Ce^{-t}}$ . For $C > 0$ these rise from $0$ to $1$ ; together with $x \equiv 0$ and $x \equiv 1$ they are all the solutions in $[0, 1]$ .
Blow-up $x' = x^2$ , $x(0) = 1$ . Here $-\frac1x = t + c$ , so $x(t) = \dfrac{1}{1 - t}$ , which tends to $+\infty$ as $t \to 1^-$ . The right-hand side $x^2$ is perfectly smooth, yet the solution exists only on $(-\infty, 1)$ .
Non-uniqueness $x' = \sqrt{|x|}$ , $x(0) = 0$ . Besides $x \equiv 0$ (from the zero of $h$ ), separation on $x > 0$ gives $2\sqrt{x} = t + c$ , and gluing gives the solution $x(t) = t^2/4$ for $t \ge 0$ , $x(t) = 0$ for $t < 0$ . Indeed for every $c \ge 0$ , the function equal to $0$ for $t \le c$ and $(t - c)^2/4$ for $t \ge c$ is a solution: infinitely many solutions share the same initial value.

x′ = √|x|: many solutions through one point

The parabola through (2, 1) reaches the axis at t = 0 — exactly where the zero solution also passes. Which way should a solution starting at (0, 0) go?

Examples

dy/dx = sqrt(abs(y))

Click the field to draw a solution, or type a start point below.

Start pointx₀y₀

Solution curves: 3. They are numerical approximations (fourth-order Runge–Kutta), not exact formulas.

Latest curve through (2, 1).

Theorem— Theorem (first-order linear equations)

Let $p, q$ be continuous on an interval $I$ , $t_0 \in I$ , and $P(t) = \int_{t_0}^t p(s)\,ds$ . The initial value problem $y' + p(t)\,y = q(t)$ , $y(t_0) = y_0$ has exactly one solution on $I$ , namely

y(t) = e^{-P(t)}\left(y_0 + \int_{t_0}^{t} q(s)\,e^{P(s)}\,ds\right).

Proof

The key identity, from the product and chain rules and $P' = p$ , is

\big(y\,e^{P}\big)' = (y' + p\,y)\,e^{P} .

So $y$ solves the equation if and only if $(y e^P)' = q e^P$ . Existence: the displayed $y$ has $y e^{P} = y_0 + \int_{t_0}^t q e^{P}$ , whose derivative is $q e^P$ by the fundamental theorem of calculus, and $y(t_0) = y_0$ . Uniqueness: if $y_1, y_2$ are solutions, $z = y_1 - y_2$ satisfies $(z e^{P})' = 0$ on the interval $I$ , so $z e^P$ is constant, equal to its value $0$ at $t_0$ .

∎

The factor $e^{P}$ is the integrating factor. With $p = -a$ and $q = 0$ the theorem says that the only solutions of $x' = ax$ are $x(t) = x(0)e^{at}$ — the law of exponential growth and decay.

Linear equations with constant coefficients. For $y'' + b\,y' + c\,y = 0$ , try $y = e^{rt}$ : it is a solution exactly when $r^2 + br + c = 0$ , the characteristic equation.

Theorem— Theorem (second-order constant-coefficient equations)

The real solutions of $y'' + by' + cy = 0$ form a two-dimensional vector space, with basis

$e^{r_1 t},\; e^{r_2 t}$ if the characteristic roots $r_1 \ne r_2$ are real;
$e^{rt},\; t\,e^{rt}$ if there is a double root $r$ ;
$e^{\alpha t}\cos\beta t,\; e^{\alpha t}\sin\beta t$ if the roots are $\alpha \pm i\beta$ with $\beta \ne 0$ .

Proof

Write $D = d/dt$ ; the equation is $(D - r_1)(D - r_2)\,y = 0$ , where $r_1, r_2$ are the roots (the factors commute because the coefficients are constant). Each listed function is a solution: $(D - r)e^{rt} = 0$ ; for a double root, $(D - r)(t e^{rt}) = e^{rt}$ and then $(D - r)e^{rt} = 0$ ; for complex roots, $e^{(\alpha + i\beta)t}$ is a complex solution, and since the equation has real coefficients its real and imaginary parts are real solutions. They are independent: their Wronskians at $t = 0$ are $r_2 - r_1$ , $1$ and $\beta$ , all non-zero. The space is two-dimensional: by the existence and uniqueness theorem below, the linear map $y \mapsto (y(0), y'(0))$ from solutions to $\R^2$ is a bijection.

∎

The damped oscillator $m y'' + \gamma y' + k y = 0$ ( $m, k > 0$ , $\gamma \ge 0$ ) has roots $\frac{-\gamma \pm \sqrt{\gamma^2 - 4mk}}{2m}$ : overdamped if $\gamma^2 > 4mk$ (two decaying exponentials), critically damped if $\gamma^2 = 4mk$ , and underdamped if $\gamma^2 < 4mk$ (decaying oscillations $e^{-\gamma t/2m}\cos(\omega t + \phi)$ ). For an inhomogeneous equation $y'' + by' + cy = g(t)$ , the general solution is one particular solution plus the general homogeneous solution, because the difference of two solutions solves the homogeneous equation.

ProofH

Linear systems and phase portraits

For a constant $d \times d$ matrix $A$ , define the matrix exponential

e^{tA} = \sum_{k=0}^{\infty} \frac{t^k A^k}{k!} = I + tA + \frac{t^2A^2}{2} + \cdots .

The series converges absolutely for every $t$ in any operator norm, since $\|A^k\| \le \|A\|^k$ , and it can be differentiated term by term (it is a power series in $t$ with infinite radius of convergence): $\frac{d}{dt}e^{tA} = A\,e^{tA} = e^{tA}A$ . If $AB = BA$ then $e^{A + B} = e^A e^B$ (multiply the series, exactly as for numbers); in particular $e^{-tA}e^{tA} = I$ .

Theorem— Theorem (linear systems)

For every $x_0 \in \R^d$ , the unique solution of $x' = Ax$ with $x(0) = x_0$ is $x(t) = e^{tA}x_0$ , defined for all $t \in \R$ . If $A$ has a basis of eigenvectors $v_1, \dots, v_d$ with eigenvalues $\lambda_1, \dots, \lambda_d$ and $x_0 = \sum c_i v_i$ , then $x(t) = \sum c_i e^{\lambda_i t} v_i$ .

Proof

$x(t) = e^{tA}x_0$ satisfies $x' = Ae^{tA}x_0 = Ax$ and $x(0) = x_0$ . If $y$ is any solution, then $\big(e^{-tA}y(t)\big)' = -Ae^{-tA}y + e^{-tA}Ay = 0$ (as $A$ commutes with $e^{-tA}$ ), so $e^{-tA}y(t) = y(0) = x_0$ and $y(t) = e^{tA}x_0$ . Finally, $A^k v_i = \lambda_i^k v_i$ gives $e^{tA}v_i = e^{\lambda_i t}v_i$ .

∎

Complex eigenvalues $\alpha \pm i\beta$ produce rotation combined with growth or decay at rate $e^{\alpha t}$ : for $A = \begin{pmatrix} 0 & 1 \\ -1 & 0 \end{pmatrix}$ (that is, $x'' = -x$ ), $e^{tA} = \begin{pmatrix} \cos t & \sin t \\ -\sin t & \cos t \end{pmatrix}$ , a rotation. Non-diagonalisable matrices, via the Jordan form, add polynomial factors $t^k e^{\lambda t}$ .

Theorem— Theorem (classification of planar linear systems)

Let $A$ be a real $2 \times 2$ matrix with trace $\tau$ and determinant $\Delta \ne 0$ , and consider $x' = Ax$ , whose only equilibrium is $0$ . Then:

$\Delta < 0$ : a saddle (unstable);
$\Delta > 0$ and $\tau^2 > 4\Delta$ : a node, asymptotically stable if $\tau < 0$ , unstable if $\tau > 0$ ;
$\Delta > 0$ and $\tau^2 < 4\Delta$ : a spiral (focus), asymptotically stable if $\tau < 0$ , unstable if $\tau > 0$ , and a centre (stable, not asymptotically) if $\tau = 0$ .

Proof

The eigenvalues are the roots of $\lambda^2 - \tau\lambda + \Delta = 0$ , so $\lambda_1\lambda_2 = \Delta$ and $\lambda_1 + \lambda_2 = \tau$ . If $\Delta < 0$ they are real with opposite signs: along one eigenvector solutions decay, along the other they grow. If $\Delta > 0$ and $\tau^2 > 4\Delta$ they are real, distinct, and of the same sign as $\tau$ : all solutions $c_1e^{\lambda_1 t}v_1 + c_2e^{\lambda_2 t}v_2$ decay (or all grow). If $\tau^2 < 4\Delta$ they are $\frac\tau2 \pm i\beta$ , and solutions are rotations multiplied by $e^{\tau t/2}$ . (The boundary case $\tau^2 = 4\Delta$ gives degenerate or star nodes, with the stability of the sign of $\tau$ .)

∎

A damped oscillator: x″ + 0.6x′ + 2x = 0

Trace −0.6, determinant 2: a stable spiral. Edit the equations to make a saddle, a centre or a node.

Examples

x′ = y, y′ = −2·x − 0.6·y

darker, longer arrows = faster forward in time backward in timeequilibrium (hollow = stable or centre)

Click the plane to start a trajectory, or type a start point below.

Start pointx₀y₀

Trajectories: 3. They are numerical approximations (fourth-order Runge–Kutta), not exact formulas.

Linear systemA = [01−2−0.6]

T = −0.6, D = 2, T² − 4D = −7.64; λ = −0.3 ± 1.382i

Equilibrium at (0, 0): spiral sink (stable)

Trace T and determinant D: below the T-axis, saddles; above the dashed parabola D = T²/4, spirals; between them, nodes; on the upper D-axis, centres. Left of the D-axis is stable, right is unstable.

ProofH

Existence and uniqueness: the Picard–Lindelöf theorem

We now prove that well-behaved initial value problems have exactly one solution, at least for a short time. We write $|\cdot|$ for a norm on $\R^d$ .

Theorem— Theorem (Picard–Lindelöf)

Let $f(t, x)$ be continuous on the box $R = \{(t, x) : |t - t_0| \le a,\ |x - x_0| \le b\}$ , with $|f| \le M$ on $R$ , and Lipschitz in $x$ on $R$ with constant $K$ . Let $h = \min(a, b/M)$ . Then the initial value problem $x' = f(t, x)$ , $x(t_0) = x_0$ has exactly one solution on $[t_0 - h, t_0 + h]$ .

Proof

Write $I = [t_0 - h, t_0 + h]$ .

Step 1: an integral equation. By the fundamental theorem of calculus, a continuous $x : I \to \R^d$ with $|x(t) - x_0| \le b$ solves the IVP if and only if

x(t) = x_0 + \int_{t_0}^{t} f(s, x(s))\,ds \qquad \text{for all } t \in I.

Step 2: Picard iterates. Put $x_0(t) = x_0$ and $x_{k+1}(t) = x_0 + \int_{t_0}^t f(s, x_k(s))\,ds$ . By induction each $x_k$ is continuous on $I$ with $|x_k(t) - x_0| \le M|t - t_0| \le Mh \le b$ , so $(s, x_k(s))$ stays in $R$ and the next iterate is defined.

Step 3: an estimate. We claim $|x_{k+1}(t) - x_k(t)| \le \dfrac{M K^k |t - t_0|^{k+1}}{(k + 1)!}$ . For $k = 0$ this is $|x_1(t) - x_0| \le M|t - t_0|$ . If it holds for $k$ , then by the Lipschitz condition (for $t \ge t_0$ ; $t \le t_0$ is symmetric)

|x_{k+2}(t) - x_{k+1}(t)| \le \int_{t_0}^{t} K\,|x_{k+1}(s) - x_k(s)|\,ds \le \int_{t_0}^t \frac{M K^{k+1}(s - t_0)^{k+1}}{(k+1)!}\,ds = \frac{M K^{k+1}(t - t_0)^{k+2}}{(k+2)!}.

Step 4: convergence. Since $\sum_k M K^k h^{k+1}/(k+1)!$ converges (to $\frac{M}{K}(e^{Kh} - 1)$ , or $Mh$ if $K = 0$ ), the Weierstrass M-test shows that $x_k = x_0 + \sum_{j < k}(x_{j+1} - x_j)$ converges uniformly on $I$ to a continuous function $x$ with $|x(t) - x_0| \le b$ . As $|f(s, x_k(s)) - f(s, x(s))| \le K|x_k(s) - x(s)| \to 0$ uniformly, we may pass to the limit under the integral sign in the definition of $x_{k+1}$ : $x$ satisfies the integral equation, hence solves the IVP.

Step 5: uniqueness (Grönwall's argument). A solution $y$ on $I$ stays in the box: up to the first time it could leave, $|y(t) - x_0| \le M|t - t_0| \le Mh \le b$ , so it never does. So let $x, y$ be two solutions with values in the box, and for $t \in [t_0, t_0 + h]$ put $u(t) = |x(t) - y(t)|$ . Subtracting the integral equations, $u(t) \le K\int_{t_0}^t u(s)\,ds$ . Let $U(t) = \int_{t_0}^t u$ ; then $U' = u \le KU$ , so

\big(e^{-K(t - t_0)}U(t)\big)' = e^{-K(t-t_0)}\big(u(t) - KU(t)\big) \le 0 .

Since $U(t_0) = 0$ , we get $e^{-K(t-t_0)}U(t) \le 0$ ; but $U \ge 0$ , so $U \equiv 0$ and then $u \le KU = 0$ . The same argument works for $t \le t_0$ .

∎

The same Grönwall argument gives continuous dependence on initial data: two solutions with $x(t_0) = x_0$ and $y(t_0) = y_0$ satisfy $|x(t) - y(t)| \le |x_0 - y_0|\,e^{K|t - t_0|}$ . The factor $e^{K|t - t_0|}$ can be enormous, which is the seed of chaos in nonlinear systems.

What the theorem does not say.

It is local. The solution of $x' = x^2$ , $x(0) = 1$ exists only for $t < 1$ , although $x^2$ is smooth everywhere. In general, if $f$ is locally Lipschitz on $\R \times \R^d$ , every solution extends to a maximal open interval $(\alpha, \omega)$ , and if $\omega < \infty$ then $|x(t)| \to \infty$ as $t \to \omega$ : solutions die only by escaping to infinity. Linear systems never blow up.
Continuity alone gives existence but not uniqueness. Peano's theorem (1890) says that if $f$ is merely continuous, a solution exists for a short time; the proof uses compactness (Arzelà–Ascoli) instead of a contraction. For $x' = \sqrt{|x|}$ the right-hand side is continuous but not Lipschitz near $x = 0$ , and uniqueness fails, as we saw.
It gives no formula. Liouville showed in 1841 that the Riccati equation $x' = t^2 + x^2$ cannot be solved by quadratures (elementary functions and their integrals); its solutions need Bessel functions. Most ODEs, including the pendulum $\theta'' = -\sin\theta$ and the three-body problem, have no elementary closed-form solution — hence the importance of qualitative and numerical methods.

ExplanationH

Qualitative behaviour

For an autonomous scalar equation $x' = f(x)$ with $f$ locally Lipschitz, uniqueness alone gives a complete picture. Solutions cannot cross the equilibria (zeros of $f$ ), so a solution starting between two consecutive equilibria stays between them, and since $x' = f(x)$ has a fixed sign there, it is monotone and tends to one of them (or to $\pm\infty$ ). Hence an equilibrium $x^*$ with $f'(x^*) < 0$ is asymptotically stable (arrows point towards it on both sides) and one with $f'(x^*) > 0$ is unstable. For the logistic equation $f(x) = x(1 - x)$ , $f'(0) = 1$ and $f'(1) = -1$ : the population moves away from extinction and settles at the carrying capacity.

In the plane, an equilibrium $x^*$ of $x' = f(x)$ is studied through the linearisation $y' = Df(x^*)\,y$ . The Hartman–Grobman theorem (Grobman 1959, Hartman 1960) says that if $Df(x^*)$ has no eigenvalue with zero real part (a hyperbolic equilibrium), then near $x^*$ the phase portrait is a continuous deformation of the linear one: saddles stay saddles and sinks stay sinks. A linear centre is not hyperbolic, and the nonlinear terms decide.

Another tool is a Lyapunov function: a function $V$ with a strict minimum at $x^*$ that does not increase along solutions. For the damped pendulum $\theta'' = -\sin\theta - \gamma\theta'$ , the energy $E = \tfrac12\theta'^2 + (1 - \cos\theta)$ satisfies $E' = \theta'\theta'' + \sin\theta\,\theta' = -\gamma\theta'^2 \le 0$ , so the bottom equilibrium is stable — without solving the equation. Energy is also a powerful uniqueness tool: for $x'' = -x$ , $E = x^2 + x'^2$ is constant, so a solution with $x(0) = x'(0) = 0$ stays at $0$ (this is formalised below).

Formal mathematicsF

Differential equations in Lean and Mathlib

In Mathlib, " $x$ solves $x' = f(t, x)$ " is written ∀ t, HasDerivAt x (f t (x t)) t. The lesson's module verifies explicit solutions, proves the integrating-factor uniqueness argument and the energy argument, and makes the failure modes precise.

Exponential growth: the solution, and its uniqueness by the integrating factor $e^{-at}$ — a quantity with zero derivative is constant:

exp_solvesLean 4Proved here from basic factsFormally verified

/-- `x(t) = x₀ e^{a t}` solves `x' = a x` with `x(0) = x₀`. -/
theorem exp_solves (a x₀ t : ℝ) :
    HasDerivAt (fun s => x₀ * exp (a * s)) (a * (x₀ * exp (a * t))) t := by
  have h := ((hasDerivAt_id' t).const_mul a).exp.const_mul x₀
  convert h using 1
  ring