The Courant-Friedrichs-Lewy condition

Courant, Friedrichs and Lewy (1928) observed that an explicit scheme for \(u_t + c\,u_x = 0\) computes each new value from a few neighbors, so information can travel at most one cell per step. If the true wave moves farther than that – Courant number \(\nu = c\,\Delta t/\Delta x > 1\) – the scheme cannot possibly see where the solution came from, and it diverges. This script advects a pulse with the upwind scheme below and above \(\nu = 1\), and with Lax-Wendroff, whose limit is the same.

import matplotlib.pyplot as plt
import numpy as np

from mathematicskit.pde import AdvectionEquation1D, check_cfl

adv = AdvectionEquation1D(lambda x: np.exp(-200 * (x - 0.5) ** 2), n=200, c=1.0)

Below and above nu = 1

fig, axes = plt.subplots(1, 2, figsize=(11, 4), sharey=True)
for ax, scheme in zip(axes, ("upwind", "lax_wendroff"), strict=False):
    ax.plot(adv.x, adv.exact(1.0), "k--", lw=2, label="exact (one period)")
    for nu, color in ((0.5, "tab:blue"), (0.95, "tab:green"), (1.2, "tab:red")):
        sol = adv.solve(1.0, dt=nu * adv.dx, scheme=scheme)
        check = sol.extra["stability"]
        print(f"{scheme:>12}  nu = {nu:4.2f}  stable = {check.stable!s:5}  max |u| = {np.max(np.abs(sol.final)):.3g}")
        ax.plot(sol.x, np.clip(sol.final, -0.5, 1.5), color=color, label=f"nu = {nu}")
    ax.set_title(scheme)
    ax.set_xlabel("$x$")
    ax.set_ylim(-0.5, 1.5)
axes[0].legend(fontsize=8)
fig.suptitle("CFL: the Courant number must not exceed 1")
fig.tight_layout()

      upwind  nu = 0.50  stable = True   max |u| = 0.707
      upwind  nu = 0.95  stable = True   max |u| = 0.952
      upwind  nu = 1.20  stable = False  max |u| = 7.95e+07
lax_wendroff  nu = 0.50  stable = True   max |u| = 0.995
lax_wendroff  nu = 0.95  stable = True   max |u| = 0.999
lax_wendroff  nu = 1.20  stable = False  max |u| = 5.95e+28

The domain-of-dependence argument, in numbers

for dt in (0.004, 0.005, 0.006):
    r = check_cfl(adv.c, dt, adv.dx, "upwind")
    print(f"dt = {dt}: Courant number {r.number:.2f} (limit {r.limit:.0f}) -> {'stable' if r.stable else 'UNSTABLE'}")

plt.show()

dt = 0.004: Courant number 0.80 (limit 1) -> stable
dt = 0.005: Courant number 1.00 (limit 1) -> stable
dt = 0.006: Courant number 1.20 (limit 1) -> UNSTABLE