3 Practical Stability

We have assumed that the largest growth rate modulus must be bounded by 1 for stability; it is possible to conclude that there will be stability in the sense that solutions remain bounded, provided the largest eigenvalue of ${\rm\bf B}_1^{-1}{\rm\bf B}_0$ satisfies

$$\vert\lambda\vert\le 1+C\Delta t,\ {\rm as\ }\Delta t\to 0,$$ (117)

for some constant $C>0$, since then

$$\vert\lambda^n\vert\le (1+C\Delta t)^n\le {\rm e}^{Cn\Delta t},$$ (118)

and so

$$\vert\lambda^n\vert\le{\rm e}^{Ct}\ {\rm as}\ n\Delta t\to t \ {\rm with}\ \Delta t\to0,\ n\to\infty.$$ (119)

This is called Lax-Richtmyer stability.

However for practical purposes, if the modulus of $\lambda$ exceeds 1 for finite $\Delta t$ then the calculations invariably become too large to handle. Hence it is safer to require practical or strict stability, $\vert\lambda\vert\le 1$.

Example: Convection Diffusion

Consider an explicit discretisation of a convection-diffusion equation,

$${{\partial u}\over{\partial t}}+a{{\partial u}\over{\partial x}}=D{{\partial^2 u}\over{\partial x^2}} ,$$ (120)

$$U^{n+1}_j=U^n_j-{1\over 2}\nu (U^n_{j+1}-U^n_{j-1}) +\mu (U^n_{j+1}-2U^n_j+U^n_{j-1}),$$ (121)

where $\nu=a\Delta t/h$, $\mu=D\Delta t/h^2$. Using a Fourier mode

$$\lambda = 1-i\nu\sin kh -4\mu\sin {1\over 2}kh,$$ (122)

so if $p=\sin^2{1\over 2}kh$,

$$\vert\lambda(p)\vert^2=(1-4\mu p)^2+4\nu^2 p(1-p).$$ (123)

(a) Lax-Richtmyer Stability

If we write $\nu^2={{a^2\mu}\over D}\Delta t$, then

$$\vert\lambda\vert^2=(1-4\mu p)^2+{{4a^2\mu p(1-p)}\over D}\Delta t,$$ (124)

and provided $0\le \mu\le{1\over 2}$, we will have

$$\vert\lambda\vert^2\le 1+{{a^2\mu}\over D}\Delta t,$$ (125)

so that the scheme will be L-R stable for $\Delta t\to 0$ with $\mu=D\Delta t/h^2$ fixed. However in practice this is not sufficiently strong, for instance with $\mu =0.25$, $\nu = 0.75$, $\vert\lambda\vert^2=1+{{17}\over 4}p-{{13}\over 4}p^2$ so that growth factors of greater than 1 occur for a range of small values of $p$ leading to breakdown of a calculation.

(b) Practical Stability

If we write

$$\vert\lambda\vert^2=1-4(2\mu -\nu^2)p+4(4\mu^2 - \nu^2)p^2,$$ (126)

then $\vert\lambda(0)\vert^2=1$ and $\vert\lambda (1)\vert^2=(1-4\mu )\le 1$ provided $0\le \mu\le{1\over 2}$. Since we are just dealing with a quadratic all we need have is the slope (with respect to $p$) to negative at $p=0$ in order for $\vert\lambda\vert$ to be less than 1 over the range $0\le p\le 1$. Hence we also need $2\mu-\nu^2\ge 0$. This can be rearranged as

$${{a^\Delta t^2}\over{2h^2}}\le {{D\Delta t}\over h^2}\le {1\over 2},$$ (127)

or

$${{a\Delta t}\over h}\le \min ({2\over{Pe}},{{Pe}\over 2}),$$ (128)

where $Pe$ is the mesh Peclet number, $Pe=ah/D$.