Note
Go to the end to download the full example code or to run this example in your browser via JupyterLite.
The Courant-Isaacson-Rees upwind scheme#
Courant, Isaacson and Rees (1952) differenced \(u_x\) on the side the information comes from: for \(c > 0\), \(u_j^{n+1} = u_j^n - \nu(u_j^n - u_{j-1}^n)\). The scheme is stable for \(\nu \le 1\), exact at \(\nu = 1\), and otherwise first-order accurate, with numerical diffusion \(\tfrac12 c\,\Delta x(1-\nu)\,u_{xx}\) that smears the pulse. Differencing on the wrong side is unstable at every step size.
import matplotlib.pyplot as plt
import numpy as np
from mathematicskit.pde import AdvectionEquation1D
pulse = lambda x: np.exp(-200 * (x - 0.5) ** 2)
Numerical diffusion shrinks as nu -> 1#
adv = AdvectionEquation1D(pulse, n=200, c=1.0)
fig, ax = plt.subplots()
ax.plot(adv.x, adv.exact(1.0), "k--", lw=2, label="exact")
for nu in (0.25, 0.5, 0.9, 1.0):
sol = adv.solve(1.0, dt=nu * adv.dx, scheme="upwind")
print(f"nu = {nu:4.2f}: peak {sol.final.max():.3f}, max error {np.max(np.abs(sol.final - adv.exact(1.0))):.2e}")
ax.plot(sol.x, sol.final, label=f"upwind, nu = {nu}")
ax.set_xlabel("$x$")
ax.set_title("Upwind after one period: diffusion (1 - nu)")
ax.legend()

nu = 0.25: peak 0.632, max error 3.68e-01
nu = 0.50: peak 0.707, max error 2.93e-01
nu = 0.90: peak 0.911, max error 8.95e-02
nu = 1.00: peak 1.000, max error 6.66e-16
<matplotlib.legend.Legend object at 0x7fd7bffc1f70>
The upwind side depends on the sign of c#
left = AdvectionEquation1D(pulse, n=200, c=-1.0)
sol = left.solve(0.25, dt=0.8 * left.dx, scheme="upwind")
print(f"c = -1: max error after t = 0.25: {np.max(np.abs(sol.final - left.exact(0.25))):.2e}")
c = -1: max error after t = 0.25: 4.79e-02
First-order convergence#
print(" n error ratio")
previous = None
for n in (100, 200, 400, 800):
a = AdvectionEquation1D(lambda x: np.sin(2 * np.pi * x), n=n)
error = np.max(np.abs(a.solve(1.0, dt=0.5 * a.dx, scheme="upwind").final - a.exact(1.0)))
print(f"{n:>4} {error:.2e} {previous / error if previous else float('nan'):5.2f}")
previous = error
plt.show()
n error ratio
100 9.40e-02 nan
200 4.82e-02 1.95
400 2.44e-02 1.98
800 1.23e-02 1.99
Total running time of the script: (0 minutes 0.184 seconds)