Von Neumann stability analysis and the CFL condition

Intuition Beginner

A finite-difference scheme marches a grid of numbers forward in time. The danger is that small errors — round-off, or just the gap between the grid answer and the true one — might grow from step to step until they swamp the answer. A scheme that keeps every error bounded is called stable; one that lets some error grow without limit is unstable. The whole point of this unit is a clean test for which schemes are which.

The test rests on one idea borrowed from sound. Any pattern on the grid can be built from pure waves of different wavelengths, just as any musical chord is a sum of pure tones. A constant-coefficient scheme treats each pure wave on its own: it does not mix one wavelength into another. So you can feed in a single wave and watch what one step does to it. Each step multiplies that wave by a number, the same number every step, called the amplification factor. If that number has size at most one, the wave stays bounded; if its size is bigger than one, the wave grows each step and the scheme blows up.

The rule is then short. Run through every wavelength the grid can hold, find the amplification factor for each, and check that none of them has size above one (give or take a tiny allowance for solutions that are genuinely supposed to grow). If all pass, the scheme is stable. This is von Neumann analysis, named for the mathematician who turned the wave-by-wave idea into a working tool.

There is a second, simpler limit that applies to waves and signals that travel — like a pulse moving down a string from the earlier wave-equation unit 02.13.04. Information in the true problem moves at a fixed speed, so in one time step it travels a fixed distance. A marching scheme also gathers information from a fixed spread of neighbours each step. If the true signal moves farther in one step than the scheme can reach, the scheme is trying to predict something it never looked at, and no choice of arithmetic can save it. The requirement that the scheme's reach cover the true signal's travel is the Courant-Friedrichs-Lewy condition, or CFL for short.

Visual Beginner

There are two pictures here. The first is the wave-by-wave test: feed one pure wave into one step and see how its height changes. The second is the travel-distance test: compare how far the true signal moves with how far the scheme can reach.

The table collects the amplification factors of the standard schemes from the parabolic unit 43.11.02, written in terms of the single number $r$ (time step over squared spacing) and the wave. The size of that factor is what decides stability.

scheme	what one step multiplies a wave by	stays at size $\le 1$ when
forward Euler (explicit heat)	$1$ minus a positive amount that can exceed $2$	only if $r\le 1/2$
backward Euler (implicit heat)	$1$ over $1$ plus a positive amount	always
Crank-Nicolson (implicit heat)	a ratio whose top is smaller than its bottom	always

 time
  ^   physical travel in one step          scheme's reach in one step
  |     (true signal moves a*k)              (stencil sees +/- one cell)
  |
  |        \                                       o
  |         \   characteristic                    /|\
  |          \  slope 1/a                         / | \
  |           *  point we predict                o  o  o   known cells
  |
  +------------------------------------------------------------> space
   CFL ok:  a*k  <=  h   (true travel fits inside the scheme's reach)
   CFL bad: a*k  >   h   (signal comes from a cell the scheme never read)

The takeaway: von Neumann analysis is the wave-by-wave size check on the amplification factor, and it gives sharp stability limits like $r\le 1/2$ for the explicit heat scheme. The CFL condition is the coarser travel-distance check; it must hold for any explicit hyperbolic scheme, but on its own it does not promise stability.

Worked example Beginner

Let us run the wave-by-wave test on the explicit heat scheme and watch it produce the limit $r\le 1/2$ with a real number. The scheme updates each grid value by $U_j^{n+1}=U_j^n+r(U_{j-1}^n-2U_j^n+U_{j+1}^n)$, where $r=k/h^2$.

Feed in a single wave. Take the most demanding one the grid can carry: a sawtooth that flips sign at every grid point, $U_j^n=(-1{)}^j$. Its neighbours are the negatives of its own value: if $U_j^n=+1$, both $U_{j-1}^n$ and $U_{j+1}^n$ equal $-1$.

Step 1. Compute the neighbour comparison at a point where $U_j^n=1$: the left neighbour plus the right neighbour minus twice the centre is $(-1)+(-1)-2(1)=-4$.

Step 2. Apply the update: $U_j^{n+1}=1+r(-4)=1-4r$. So one step multiplied this wave by the factor $g=1-4r$. Every point on the sawtooth gets multiplied by the same $g$, which is the whole point of the wave-by-wave view.

Step 3. Ask when the size of $g$ stays at most one: $|1-4r|\le 1$. This means $-1\le 1-4r\le 1$, that is $0\le 4r\le 2$, so $r\le 1/2$.

Step 4. See the failure past the limit. Take $r=1$, twice the limit. Then $g=1-4=-3$, of size three. The sawtooth flips sign and triples in height every step: $1,-3,9,-27,\dots$. The error explodes.

What this tells us: the single worst wave, the sign-flipping sawtooth, already pins down the exact stability limit $r\le 1/2$ — the same limit the parabolic unit 43.11.02 got from eigenvalues, now from a one-line wave test. The sawtooth is worst because the neighbour comparison hits it hardest, which is why von Neumann analysis always checks every wavelength and reports the most dangerous one.

Check your understanding Beginner

Formal definition Intermediate+

Work on a uniform spatial grid $x_j=jh$, $j\in \mathbb{Z}$, with time step $k$ and time levels $t_n=nk$, and consider a constant-coefficient linear one-step scheme on a periodic domain (or the whole line), $$ U^{n+1}j = \sum{\ell = -p}^{q} c_\ell, U^n_{j+\ell}, $$ with fixed real (or complex) coefficients $c_{\ell }$ independent of $j$ and $n$; the parabolic FTCS, backward-Euler, and Crank-Nicolson schemes of 43.11.02 and the advection schemes of the next unit are all of this form (implicit schemes after solving the constant-coefficient recurrence). Because the coefficients do not depend on $j$, the scheme commutes with the grid shift, and the discrete Fourier transform diagonalises it.

The Fourier mode and the amplification factor. Insert the single Fourier mode $$ U^n_j = g^n, e^{i\xi j h}, \qquad \xi \in [-\pi/h, \pi/h], $$ where $\xi$ is the spatial wavenumber and $g$ is to be determined. Substituting and dividing by the common factor $g^n{e}^{i\xi jh}$ leaves a single scalar equation for $g$. For an explicit scheme this gives $g$ outright; for an implicit scheme it gives $g$ after dividing the contributions from the new level. The resulting $$ g(\xi) = \sum_{\ell=-p}^{q} c_\ell, e^{i\xi \ell h} $$ (explicit case) is the amplification factor — the symbol of the scheme, the complex number by which the scheme multiplies the mode of wavenumber $\xi$ in one step. After $n$ steps that mode is multiplied by $g(\xi {)}^n$.

The von Neumann stability condition. The scheme is von Neumann stable (equivalently, ${\ell }^2$-stable in the Lax-Richtmyer sense of 43.11.05) if there is a constant $C$, independent of $\xi$, $h$, $k$, such that $$ |g(\xi)| \le 1 + C k \qquad \text{for all } \xi \in [-\pi/h, \pi/h]. $$ When the exact solution does not grow (the parabolic and pure-advection cases), the sharp form is the strict condition $|g(\xi )|\le 1$ for all $\xi$. The $Ck$ slack permits genuine exponential growth at rate $O(1)$ in the true solution: over $n\le T/k$ steps, $(1+Ck{)}^n\le e^{CT}$ stays bounded. For a system of $s$ equations the coefficients $c_{\ell }$ are $s\times s$ matrices, $g$ becomes the amplification matrix $G(\xi )\in {\mathbb{C}}^{s\times s}$, and the condition is uniform power-boundedness $\Vert G(\xi {)}^n\Vert \le K$ for $nk\le T$ — the scalar $|g|\le 1+Ck$ is necessary but, for non-normal $G$, not sufficient.

The Courant-Friedrichs-Lewy (CFL) condition. For a hyperbolic problem with finite propagation speed, the analytic domain of dependence of a point $(x,t)$ is the set of initial points whose data determine $u(x,t)$; for the advection equation $u_t+a{u}_x=0$ of 02.13.04 it is the single foot of the characteristic, $x-at$. The numerical domain of dependence is the set of initial grid points that the scheme's stencil, iterated to time $t$, actually consults. The CFL condition is the requirement that the numerical domain of dependence contain the analytic one. For an explicit scheme with stencil width one cell per step, applied to $u_t+a{u}_x=0$, this is $$ \nu := \frac{|a|,k}{h} \le 1, $$ where $\nu$ is the Courant number. The symbols $g(\xi )$, $\xi$, the amplification matrix $G(\xi )$, the Courant number $\nu$, the mesh ratio $r=k/h^2$, and the wavenumber range $[-\pi /h,\pi /h]$ are recorded in _meta/NOTATION.md.

Counterexamples to common slips

• The CFL condition is necessary, not sufficient. It is a statement about which grid points the scheme reads, derived from geometry alone; it knows nothing about how the scheme weights them. The centred explicit scheme for advection satisfies CFL for $\nu \le 1$ yet is von Neumann unstable for every $\nu >0$. CFL screens out hopeless schemes; von Neumann analysis decides the survivors.

• Von Neumann analysis assumes constant coefficients and periodicity. The Fourier mode is an exact eigenfunction of the shift only when the coefficients $c_{\ell }$ are constant and the domain is periodic or the whole line. For variable coefficients one freezes the coefficients locally (a heuristic), and for non-periodic boundaries the boundary closure can destabilise a von Neumann-stable interior scheme — the failure the GKS theory addresses.

• $|g(\xi )|\le 1$ pointwise is not the same as $\Vert G(\xi {)}^n\Vert$ bounded for systems. A non-normal amplification matrix can have spectral radius $\le 1$ while its powers transiently grow; the scalar eigenvalue condition then misses an instability. The correct systems condition is uniform power-boundedness, characterised by the Kreiss matrix theorem.

• The worst mode is usually the sawtooth, not the smoothest. For diffusive schemes the binding wavenumber is $\xi h=\pi$ (the sign-flipping mode $e^{i\pi j}=(-1{)}^j$), where ${\sin }^2(\xi h/2)=1$; checking only low wavenumbers $\xi h\approx 0$, where $g\approx 1$ for any consistent scheme, certifies nothing.

Key theorem with proof Intermediate+

The signature result is that, for constant-coefficient one-step schemes on a periodic grid, von Neumann analysis is exactly ${\ell }^2$ stability: the Fourier transform turns the one-step operator into multiplication by $g(\xi )$, so the operator norm is the supremum of $|g|$, and stability is the uniform bound on that supremum. Specialised to the parabolic schemes this recovers the restriction $r\le 1/2$ and the unconditional stability of the implicit schemes by a Fourier route, matching the eigenvalue computation of 43.11.02; specialised to advection it produces the CFL number ^{[LeVeque, R. J. — Finite Difference Methods for Ordinary and Partial Differential Equations: Steady-State and Time-Dependent Problems; Strikwerda, J. C. — Finite Difference Schemes and Partial Differential Equations (2nd ed.)]}.

Theorem (von Neumann analysis equals ${\ell }^2$ stability; the heat amplification factors). Let $U_j^{n+1}={\sum }_{\ell }{c}_{\ell }{U}_{j+\ell }^n$ be a constant-coefficient one-step scheme on the periodic grid, with one-step operator $C$ and amplification factor $g(\xi )={\sum }_{\ell }{c}_{\ell }{e}^{i\xi \ell h}$. Then $\Vert C{\Vert }_{{\ell }^2}={\sup }_{\xi }|g(\xi )|$, and the scheme is ${\ell }^2$-stable iff ${\sup }_{\xi }|g(\xi )|\le 1+Ck$ for some $C$ independent of $h,k$. In particular, for $u_t=u_{xx}$ with $r=k/h^2$:

(i) FTCS has $g(\xi )=1-4r{\sin }^2(\xi h/2)$, and ${\sup }_{\xi }|g|\le 1⟺r\le 1/2$.

(ii) Backward Euler has $g(\xi )=(1+4r{\sin }^2(\xi h/2){)}^{-1}$ and Crank-Nicolson has $g(\xi )=\frac{1-2r{\sin }^2(\xi h/2)}{1+2r{\sin }^2(\xi h/2)}$; both satisfy $|g(\xi )|\le 1$ for all $\xi$ and all $r>0$.

Proof. On ${\ell }^2(\mathbb{Z})$ the discrete Fourier transform $\hat{V}(\xi )=h{\sum }_j{V}_j{e}^{-i\xi x_j}$ is, by Parseval, a unitary map onto $L^2([-\pi /h,\pi /h])$ (up to the fixed normalisation). The grid shift $(SV{)}_j=V_{j+1}$ becomes multiplication by $e^{i\xi h}$: $\widehat{SV}(\xi )=e^{i\xi h}\hat{V}(\xi )$. Since $C={\sum }_{\ell }{c}_{\ell }{S}^{\ell }$, the transform turns $C$ into multiplication by $g(\xi )={\sum }_{\ell }{c}_{\ell }{e}^{i\xi \ell h}$: $$ \widehat{CV}(\xi) = g(\xi),\widehat{V}(\xi). $$ A multiplication operator on $L^2$ has norm equal to the essential supremum of its multiplier, so $\Vert C{\Vert }_{{\ell }^2}=\Vert g{\Vert }_{L^{\infty }}={\sup }_{\xi }|g(\xi )|$, and $\Vert C^n{\Vert }_{{\ell }^2}={\sup }_{\xi }|g(\xi ){|}^n$. Uniform power-boundedness $\Vert C^n\Vert \le K$ for $nk\le T$ holds iff ${\sup }_{\xi }|g(\xi )|\le 1+Ck$: if the bound holds then ${\sup }_{\xi }|g{|}^n\le (1+Ck{)}^{T/k}\le e^{CT}$, and conversely if ${\sup }_{\xi }|g({\xi }_0)|=1+\delta$ with $\delta /k\to \infty$ then $|g({\xi }_0){|}^n$ is unbounded over $nk\le T$. This is the von Neumann condition.

(i) FTCS. The mode $U_j^n=g^n{e}^{i\xi jh}$ satisfies $g{e}^{i\xi jh}=e^{i\xi jh}+r(e^{i\xi (j-1)h}-2e^{i\xi jh}+e^{i\xi (j+1)h})$. Divide by $e^{i\xi jh}$: $$ g = 1 + r(e^{-i\xi h} - 2 + e^{i\xi h}) = 1 + 2r(\cos\xi h - 1) = 1 - 4r\sin^2(\xi h/2), $$ using $1-\cos \theta =2{\sin }^2(\theta /2)$. Since ${\sin }^2(\xi h/2)$ ranges over $[0,1]$ as $\xi$ runs over $[-\pi /h,\pi /h]$, $g$ ranges over $[1-4r,1]$. The upper end is always $\le 1$; the lower end satisfies $1-4r\ge -1⟺r\le 1/2$. Thus ${\sup }_{\xi }|g|\le 1$ exactly when $r\le 1/2$, the binding mode being $\xi h=\pi$ (the sawtooth).

(ii) Backward Euler and Crank-Nicolson. For backward Euler, $g{e}^{i\xi jh}=e^{i\xi jh}+rg(e^{-i\xi h}-2+e^{i\xi h})e^{i\xi jh}$, so $g=1+rg(2\cos \xi h-2)$, giving $g(1+4r{\sin }^2(\xi h/2))=1$ and $g=1/(1+4r{\sin }^2(\xi h/2))\in (0,1]$, so $|g|\le 1$ for all $\xi$ and all $r>0$. For Crank-Nicolson the spatial operator is split equally across the two levels; writing $s=2r{\sin }^2(\xi h/2)\ge 0$, the recurrence $g(1+s)=(1-s)$ gives $g=(1-s)/(1+s)$, and $|1-s|\le 1+s$ for every $s\ge 0$, so $|g|\le 1$ for all $\xi$, $r$. Both are unconditionally stable. $\square$

Bridge. This theorem is the foundational reason von Neumann analysis works at all: the Fourier transform is exactly the change of basis that diagonalises every constant-coefficient scheme at once, so the operator norm collapses to a supremum of scalars $|g(\xi )|$, and stability becomes a pointwise check on the symbol. This is exactly the eigen-decoupling of 43.11.02 carried from the finite sine basis on an interval to the continuous Fourier basis on the line — the discrete eigenvalues ${\mu }_p=-\frac{4}{h^2}{\sin }^2(p\pi h/2)$ there reappear here as $k\mu (\xi )=-4r{\sin }^2(\xi h/2)$ inside $g=1+k\mu (\xi )$, and the same restriction $r\le 1/2$ falls out. The Fourier-stability half is dual to the consistency half: one bounds how the scheme amplifies each mode, the other bounds the residual the exact solution leaves in that mode, and putting these together is the central insight that a consistent, von Neumann-stable scheme converges. This builds toward the CFL condition below and appears again in the Lax-Richtmyer equivalence theorem 43.11.05, where $\Vert C^n{\Vert }_{{\ell }^2}={\sup }_{\xi }|g(\xi ){|}^n$ is precisely the power-boundedness whose uniformity is the stability hypothesis, the bridge being that von Neumann analysis is the constant-coefficient engine that verifies that hypothesis mode by mode.

Exercises Intermediate+

Advanced results Master

Theorem 1 (Fourier diagonalisation and the amplification factor). On the periodic grid the discrete Fourier transform is a unitary isomorphism ${\ell }^2\to L^2([-\pi /h,\pi /h])$ carrying the grid shift to multiplication by $e^{i\xi h}$, so any constant-coefficient one-step operator $C={\sum }_{\ell }{c}_{\ell }{S}^{\ell }$ becomes multiplication by its symbol $g(\xi )={\sum }_{\ell }{c}_{\ell }{e}^{i\xi \ell h}$. Consequently $\Vert C^n{\Vert }_{{\ell }^2}={\sup }_{\xi }|g(\xi ){|}^n$, and von Neumann stability ${\sup }_{\xi }|g(\xi )|\le 1+Ck$ is identical to ${\ell }^2$ Lax-Richtmyer stability of 43.11.05. For an $s$-component system, $g$ is replaced by the amplification matrix $G(\xi )$ and the criterion is uniform power-boundedness $\Vert G(\xi {)}^n\Vert \le K$ for $nk\le T$ ^{[Richtmyer, R. D. & Morton, K. W. — Difference Methods for Initial-Value Problems (2nd ed.)]}.

Theorem 2 (the heat and advection symbols and their stability thresholds). For $u_t=u_{xx}$ with $r=k/h^2$: FTCS has $g=1-4r{\sin }^2(\xi h/2)$, stable iff $r\le 1/2$; backward Euler $g=(1+4r{\sin }^2(\xi h/2){)}^{-1}$ and Crank-Nicolson $g=\frac{1-2r{\sin }^2(\xi h/2)}{1+2r{\sin }^2(\xi h/2)}$, unconditionally stable. For $u_t+a{u}_x=0$ with $\nu =ak/h$: upwind $|g{|}^2=1-2\nu (1-\nu )(1-\cos \xi h)$, Lax-Friedrichs $|g{|}^2=1-(1-{\nu }^2){\sin }^2\xi h$, Lax-Wendroff $|g{|}^2=1-{\nu }^2(1-{\nu }^2)(1-\cos \xi h{)}^2$, all stable iff $|\nu |\le 1$, while the centred explicit scheme has $|g{|}^2=1+{\nu }^2{\sin }^2\xi h>1$ and is unconditionally unstable ^{[LeVeque, R. J. — Finite Difference Methods for Ordinary and Partial Differential Equations: Steady-State and Time-Dependent Problems]}.

Theorem 3 (CFL necessity via domains of dependence). For an explicit scheme approximating a hyperbolic problem with characteristic speed $a$, let ${\mathcal{D}}_{\mathrm{num}}(x,t)$ and ${\mathcal{D}}_{\mathrm{an}}(x,t)$ be the numerical and analytic domains of dependence. If ${\mathcal{D}}_{\mathrm{an}}(x,t)⊈{\mathcal{D}}_{\mathrm{num}}(x,t)$, then the scheme cannot converge: the initial data can be perturbed inside ${\mathcal{D}}_{\mathrm{an}}\setminus {\mathcal{D}}_{\mathrm{num}}$, changing the true value $u(x,t)$ while leaving every grid value the scheme consults unchanged. For the one-cell explicit advection stencil this is $\nu =|a|k/h\le 1$. The CFL condition is therefore a necessary convergence condition independent of the scheme's coefficients; it is not sufficient, the centred scheme being the standard witness ^{[Strikwerda, J. C. — Finite Difference Schemes and Partial Differential Equations (2nd ed.)]}.

Theorem 4 (the Kreiss matrix theorem and the systems case). For a one-step scheme for a first-order system, von Neumann's scalar condition $\rho (G(\xi ))\le 1+Ck$ (spectral radius) is necessary but, when $G(\xi )$ is non-normal, not sufficient: transient growth $\Vert G^n\Vert \gg \rho (G{)}^n$ can occur. The Kreiss matrix theorem characterises uniform power-boundedness ${\sup }_{\xi }{\sup }_n\Vert G(\xi {)}^n\Vert <\infty$ by four equivalent conditions — a resolvent bound $\Vert (zI-G{)}^{-1}\Vert \le K/(|z|-1)$ for $|z|>1$, a numerical-range condition, a Lyapunov-function condition, and a uniform-diagonalisability condition with bounded transformation — making the systems stability problem decidable in principle. The scalar von Neumann condition is the special case where $G(\xi )$ is normal (in particular for the heat schemes, where $g$ is a single eigenvalue) ^{[Richtmyer, R. D. & Morton, K. W. — Difference Methods for Initial-Value Problems (2nd ed.)]}.

Theorem 5 (GKS theory for initial-boundary-value problems). Von Neumann analysis is exact only for periodic or whole-line problems; a numerical boundary closure can destabilise a von Neumann-stable interior scheme. The Gustafsson-Kreiss-Sundström (GKS) theory replaces the Fourier mode by a normal-mode analysis in which one seeks spatially decaying solutions $U_j^n=z^n{\kappa }^j$ (with $|z|>1$, $|\kappa |<1$) of the discrete scheme satisfying the homogeneous numerical boundary conditions; the absence of such generalised eigenvalues on and outside the unit circle (the Kreiss condition) is necessary and sufficient for boundary stability. The interior von Neumann condition and the GKS boundary condition together give stability for the non-periodic initial-boundary-value problem ^{[Strikwerda, J. C. — Finite Difference Schemes and Partial Differential Equations (2nd ed.)]}.

Synthesis. The foundational reason von Neumann analysis is the working stability tool is a single structural fact: on a periodic grid the Fourier transform simultaneously diagonalises every constant-coefficient scheme, so the ${\ell }^2$ operator norm collapses to the supremum of a scalar symbol $g(\xi )$ and stability becomes a pointwise inequality. This is exactly the eigen-decoupling of 43.11.02 — the finite sine basis on an interval becoming the continuous Fourier basis on the line — and it is dual to the consistency analysis: the symbol $g(\xi )=1+k\mu (\xi )+O(k^2)$ carries both the leading consistency error in its phase and the stability budget in its modulus, so a single computation governs both halves of convergence. The central insight is that two stability statements of different strength coexist, the sharp spectral one (von Neumann, the modulus of $g$) and the coarse geometric one (CFL, the reach of the stencil), and putting these together explains the hyperbolic schemes: CFL screens which Courant numbers are even admissible, and von Neumann decides which admissible schemes survive — the centred scheme passing the first test and failing the second. The bridge to the wider theory is that von Neumann analysis is the constant-coefficient computational engine of the Lax-Richtmyer equivalence theorem 43.11.05: it verifies the power-boundedness hypothesis $\Vert C^n{\Vert }_{{\ell }^2}\le K$ mode by mode, just as the explicit spectral decomposition of 43.11.02 verifies it eigenvalue by eigenvalue, and the Kreiss matrix theorem and GKS theory generalise this verification to systems and to non-periodic boundaries.

Full proof set Master

Proposition 1 (operator norm equals the symbol supremum). Let $C={\sum }_{\ell =-p}^q{c}_{\ell }{S}^{\ell }$ on ${\ell }^2(\mathbb{Z})$, with $S$ the unit shift and $c_{\ell }\in \mathbb{C}$. Then $C$ is bounded with $\Vert C{\Vert }_{{\ell }^2}={\sup }_{\xi \in [-\pi /h,\pi /h]}|g(\xi )|$, where $g(\xi )={\sum }_{\ell }{c}_{\ell }{e}^{i\xi \ell h}$, and $\Vert C^n{\Vert }_{{\ell }^2}={\sup }_{\xi }|g(\xi ){|}^n$.

Proof. The discrete Fourier transform $\mathcal{F}:{\ell }^2(\mathbb{Z})\to L^2([-\pi /h,\pi /h])$, $(\mathcal{F}V)(\xi )=h^{1/2}{\sum }_j{V}_j{e}^{-i\xi x_j}$ (with $x_j=jh$), is unitary by Parseval's identity for Fourier series. The shift acts by $\mathcal{F}(SV)(\xi )=e^{i\xi h}(\mathcal{F}V)(\xi )$, since shifting the index by one multiplies each Fourier coefficient by $e^{i\xi h}$. Hence $\mathcal{F}C{\mathcal{F}}^{-1}$ is multiplication by $g(\xi )={\sum }_{\ell }{c}_{\ell }{e}^{i\xi \ell h}$, a bounded continuous function on the compact interval. A multiplication operator $M_g$ on $L^2$ of a $\sigma$-finite measure has $\Vert M_g\Vert =\Vert g{\Vert }_{L^{\infty }}$, and since $g$ is continuous this essential supremum is the ordinary supremum ${\sup }_{\xi }|g(\xi )|$. Unitary conjugation preserves the operator norm, so $\Vert C\Vert ={\sup }_{\xi }|g(\xi )|$. Finally $\mathcal{F}{C}^n{\mathcal{F}}^{-1}=M_{g^n}$, so $\Vert C^n\Vert ={\sup }_{\xi }|g(\xi {)}^n|={\sup }_{\xi }|g(\xi ){|}^n$. $\square$

Proposition 2 (von Neumann condition $\iff$ ${\ell }^2$ stability). The family $\{C(k)\}$ of one-step operators is ${\ell }^2$-stable — $\sup \{\Vert C(k{)}^n\Vert :0\le nk\le T,0<k\le k_0\}=:K<\infty$ — if and only if there is a constant $C^{\star }$ with ${\sup }_{\xi }|g_k(\xi )|\le 1+C^{\star }k$ for all $0<k\le k_0$.

Proof. By Proposition 1, $\Vert C(k{)}^n\Vert ={\sup }_{\xi }|g_k(\xi ){|}^n$. ($\Leftarrow$) If ${\sup }_{\xi }|g_k|\le 1+C^{\star }k$, then for $nk\le T$, $\Vert C(k{)}^n\Vert \le (1+C^{\star }k{)}^n\le (1+C^{\star }k{)}^{T/k}\le e^{C^{\star }T}=:K$, uniformly in $k$. ($\Rightarrow$) Suppose no such $C^{\star }$ exists; then there are sequences $k_m\downarrow 0$ and ${\xi }_m$ with $|g_{k_m}({\xi }_m)|=1+{\delta }_m$ and ${\delta }_m/k_m\to \infty$. Choose $n_m=\lfloor T/k_m\rfloor$, so $n_m{k}_m\to T$ and $n_m{\delta }_m\ge (T/k_m-1){\delta }_m=T{\delta }_m/k_m-{\delta }_m\to \infty$. Then $\Vert C(k_m{)}^{n_m}\Vert \ge |g_{k_m}({\xi }_m){|}^{n_m}=(1+{\delta }_m{)}^{n_m}\ge e^{n_m{\delta }_m/2}\to \infty$ (using $\log (1+\delta )\ge \delta /2$ for small $\delta$), contradicting stability. Hence the von Neumann bound is necessary. $\square$

Proposition 3 (the explicit-heat threshold is exactly $r\le 1/2$). FTCS for $u_t=u_{xx}$ has $g(\xi )=1-4r{\sin }^2(\xi h/2)$, $r=k/h^2$, and ${\sup }_{\xi }|g(\xi )|\le 1$ if and only if $r\le 1/2$.

Proof. Substituting the mode into $U_j^{n+1}=U_j^n+r(U_{j-1}^n-2U_j^n+U_{j+1}^n)$ and cancelling $g^n{e}^{i\xi jh}$ gives $g=1+r(e^{-i\xi h}-2+e^{i\xi h})=1+2r(\cos \xi h-1)=1-4r{\sin }^2(\xi h/2)$, using $\cos \theta -1=-2{\sin }^2(\theta /2)$. As $\xi$ ranges over $[-\pi /h,\pi /h]$, $w:={\sin }^2(\xi h/2)$ ranges over $[0,1]$, so $g=1-4rw\in [1-4r,1]$. The upper bound $g\le 1$ holds always (as $r,w\ge 0$). The lower bound is $1-4rw\ge -1$, i.e. $rw\le 1/2$; the worst case is $w=1$ ($\xi h=\pi$), giving $r\le 1/2$. Thus ${\sup }_{\xi }|g|=\max (1,|1-4r|)\le 1⟺4r\le 2⟺r\le 1/2$. $\square$

Proposition 4 (unconditional stability of backward Euler and Crank-Nicolson). Backward Euler has $g(\xi )=(1+4r{\sin }^2(\xi h/2){)}^{-1}$ and Crank-Nicolson has $g(\xi )=\frac{1-2r{\sin }^2(\xi h/2)}{1+2r{\sin }^2(\xi h/2)}$; both satisfy ${\sup }_{\xi }|g(\xi )|\le 1$ for every $r>0$.

Proof. Write $w={\sin }^2(\xi h/2)\in [0,1]$. Backward Euler $U_j^{n+1}-r(D^2{U}^{n+1}{)}_j=U_j^n$ gives, on the mode, $g(1+4rw)=1$, so $g=1/(1+4rw)\in (0,1]$ since $4rw\ge 0$; hence $|g|\le 1$ for all $\xi ,r$. Crank-Nicolson splits the second difference equally between levels: $g(1+2rw)=(1-2rw)$, so $g=(1-2rw)/(1+2rw)$. With $s=2rw\ge 0$, $|g|=|1-s|/(1+s)$; since $|1-s|\le 1+s$ for all $s\ge 0$ (square both sides: $(1-s{)}^2\le (1+s{)}^2⟺-2s\le 2s$, true for $s\ge 0$), $|g|\le 1$ for every $\xi$ and every $r>0$. Both stability functions are bounded by one across the whole spectrum, the Fourier counterpart of the A-stability of 43.11.02. $\square$

Proposition 5 (CFL necessity for the explicit advection stencil). Let an explicit scheme for $u_t+a{u}_x=0$ ($a>0$) compute $U_j^n$ from initial data at grid indices in $[j-n,j]$ (a one-sided one-cell stencil). If $\nu =ak/h>1$, the scheme does not converge to $u$ at fixed mesh ratio $k/h$.

Proof. The exact solution is constant along characteristics $x-at=\text{const}$, so $u(x_j,t_n)=u_0(x_j-a{t}_n)=u_0(x_j-\nu nh)$, depending only on the initial datum at the characteristic foot $\zeta :=x_j-\nu nh$. The scheme's value $U_j^n$ is a function of the initial data $\{U_i^0:j-n\le i\le j\}$, i.e. of $u_0$ on $[x_j-nh,x_j]$ — the numerical domain of dependence. If $\nu >1$ then $\zeta =x_j-\nu nh<x_j-nh$, so $\zeta$ lies strictly to the left of the numerical domain. Pick two initial data $u_0,{\tilde{u}}_0$ agreeing on $[x_j-nh,x_j]$ but differing at $\zeta$ (possible since $\zeta$ is outside that interval and $u_0$ may be chosen continuous): the scheme produces the same $U_j^n$ for both, while $u(x_j,t_n)=u_0(\zeta )\ne {\tilde{u}}_0(\zeta )=\tilde{u}(x_j,t_n)$. Refining with $k/h$ (hence $\nu$) fixed keeps $\zeta$ outside the numerical domain, so $|U_j^n-u(x_j,t_n)|$ cannot tend to zero for both data. Thus $\nu \le 1$ is necessary for convergence. $\square$

Connections Master

• The parabolic method-of-lines unit 43.11.02 supplies the test schemes whose stability this unit recomputes by Fourier modes: the FTCS, backward-Euler, and Crank-Nicolson updates are exactly the $\theta$-method schemes built there, and the amplification factor $g(\xi )=1+k\mu (\xi )$ with $k\mu (\xi )=-4r{\sin }^2(\xi h/2)$ is the continuous-spectrum version of that unit's discrete eigenvalues $k{\mu }_p=-4r{\sin }^2(p\pi h/2)$. Where 43.11.02 diagonalised the finite matrix $A_h$ by the $m$ sine modes on a bounded interval, this unit diagonalises the same constant-coefficient operator by the full continuum of Fourier modes on the line, and both recover the identical restriction $r\le 1/2$.

• The elliptic finite-difference unit 43.11.01 is the steady-state companion: the second-difference operator whose symbol is $-4r{\sin }^2(\xi h/2)/k$ here is the time-step of the very matrix $A_h$ assembled there, and the sawtooth mode $\xi h=\pi$ that binds the explicit stability limit is the discrete analogue of the highest eigenmode of $A_h$ that set the conditioning $\kappa =O(h^{-2})$. The largest-eigenvalue phenomenon that controls the parabolic time step is the same one the Fourier symbol exposes at $\xi h=\pi$.

• The Fourier-series theory 02.10.01 is the analytic engine: von Neumann analysis is the discrete Fourier transform applied to a translation-invariant operator, and the unitarity of that transform (Parseval/Plancherel) is exactly what turns the operator norm into the supremum of the symbol. The Riemann-Lebesgue decay and completeness of the Fourier basis from that unit are what justify decomposing an arbitrary grid error into modes and tracking each independently.

• The wave-equation unit 02.13.04 supplies the geometric content of the CFL condition: its characteristics and finite domains of dependence are precisely the analytic domain of dependence that the numerical domain must contain, and the d'Alembert picture of information travelling at fixed speed is what makes the travel-distance argument $ak\le h$ a hard necessary condition rather than a heuristic.

• The forthcoming hyperbolic-schemes unit 43.11.04 consumes the CFL number and the amplification-factor computations developed here: the upwind, Lax-Friedrichs, and Lax-Wendroff symbols derived in the exercises are the stability backbone of that unit, and the centred-scheme instability shown here is exactly the failure those dissipative schemes are designed to cure. The capstone Lax-Richtmyer equivalence theorem 43.11.05 is where von Neumann stability is shown to be the ${\ell }^2$ stability whose equivalence with convergence (given consistency) is the theorem's content, with the Kreiss matrix theorem providing the systems-level characterisation of the power-boundedness this unit verifies mode by mode.

Historical & philosophical context Master

The Fourier-mode stability test is named for John von Neumann, who developed it in the 1940s during wartime computational work at Los Alamos on the numerical integration of hydrodynamic and detonation problems; the method circulated through unpublished reports and was given its first widely cited exposition by John Crank, Phyllis Nicolson, and contemporaries, and in the joint analysis of John von Neumann and Robert D. Richtmyer's 1950 paper on artificial viscosity for shock computations (Journal of Applied Physics 21, 232–237) ^{[von Neumann & Richtmyer 1950]}. The systematic treatment of the amplification factor and matrix as the spectral form of stability, and its connection to the Lax-Richtmyer theory, was consolidated by Robert Richtmyer and K. W. Morton in Difference Methods for Initial-Value Problems (Interscience, 1957; 2nd ed. 1967) ^{[Richtmyer, R. D. & Morton, K. W. — Difference Methods for Initial-Value Problems (2nd ed.)]}.

The domain-of-dependence condition predates the spectral test: Richard Courant, Kurt Friedrichs, and Hans Lewy introduced it in their 1928 paper Über die partiellen Differenzengleichungen der mathematischen Physik (Mathematische Annalen 100, 32–74) ^{[Courant, Friedrichs & Lewy 1928]}, where it appeared not as a numerical recipe but as a condition for the convergence of difference approximations to hyperbolic problems, the necessity argument resting on the same characteristic-cone geometry used above. The characterisation of uniform power-boundedness that makes the systems case decidable is the Kreiss matrix theorem of Heinz-Otto Kreiss (1962), and the extension of stability analysis to non-periodic initial-boundary-value problems is the Gustafsson-Kreiss-Sundström theory of Bertil Gustafsson, Heinz-Otto Kreiss, and Arne Sundström (1972), developed in their Time-Dependent Problems and Difference Methods tradition ^{[Strikwerda, J. C. — Finite Difference Schemes and Partial Differential Equations (2nd ed.)]}.

Bibliography Master

@article{vonneumannrichtmyer1950,
  author  = {von Neumann, John and Richtmyer, Robert D.},
  title   = {A method for the numerical calculation of hydrodynamic shocks},
  journal = {Journal of Applied Physics},
  volume  = {21},
  number  = {3},
  year    = {1950},
  pages   = {232--237}
}

@article{courantfriedrichslewy1928cfl,
  author  = {Courant, Richard and Friedrichs, Kurt and Lewy, Hans},
  title   = {{\"U}ber die partiellen Differenzengleichungen der mathematischen Physik},
  journal = {Mathematische Annalen},
  volume  = {100},
  number  = {1},
  year    = {1928},
  pages   = {32--74}
}

@book{richtmyermorton1967,
  author    = {Richtmyer, Robert D. and Morton, K. W.},
  title     = {Difference Methods for Initial-Value Problems},
  edition   = {2},
  publisher = {Interscience Publishers (Wiley)},
  year      = {1967}
}

@book{strikwerda2004,
  author    = {Strikwerda, John C.},
  title     = {Finite Difference Schemes and Partial Differential Equations},
  edition   = {2},
  publisher = {Society for Industrial and Applied Mathematics (SIAM)},
  year      = {2004}
}

@book{leveque2007fdmvonneumann,
  author    = {LeVeque, Randall J.},
  title     = {Finite Difference Methods for Ordinary and Partial Differential Equations: Steady-State and Time-Dependent Problems},
  publisher = {Society for Industrial and Applied Mathematics (SIAM)},
  year      = {2007}
}

@book{gustafssonkreissoliger1995,
  author    = {Gustafsson, Bertil and Kreiss, Heinz-Otto and Oliger, Joseph},
  title     = {Time-Dependent Problems and Difference Methods},
  publisher = {John Wiley and Sons},
  year      = {1995}
}

@book{strang2007csevonneumann,
  author    = {Strang, Gilbert},
  title     = {Computational Science and Engineering},
  publisher = {Wellesley-Cambridge Press},
  year      = {2007}
}