1-D Linear Convection with the Finite Difference Method

An introduction to numerical transport, numerical diffusion, and the CFL condition using explicit upwind differencing.

The linear convection (or advection) equation is one of the most fundamental models in fluid dynamics and numerical analysis. It transports a scalar profile at a constant speed without changing its shape in the continuous form.

This simple equation is a crucial starting point because it introduces numerical transport, teaches the Courant-Friedrichs-Lewy (CFL) stability condition, and reveals important numerical errors like numerical diffusion and dispersion. It is an essential prototype to understand before tackling more complex CFD transport equations.

Applications and Context

  • Scalar transport: Modeling the movement of passive tracers (like dye in a stream) moving at a constant velocity.
  • Wave propagation prototypes: Understanding how waves travel in a medium.
  • Conservation laws: It is the advective part of the more complex Navier-Stokes conservation equations.
  • CFD method verification: It serves as a benchmark for testing new numerical schemes and teaching basic principles.

Governing Equation

The continuous 1-D linear convection equation is a partial differential equation (PDE) given by:

ut+cux=0\frac{\partial u}{\partial t} + c \frac{\partial u}{\partial x} = 0

Where:

  • uu is the scalar field we are tracking.
  • cc is the constant wave propagation speed (we assume c>0c > 0 here).
  • xx is the spatial coordinate.
  • tt is time.

Analytic Solution

For a given initial condition u(x,0)=u0(x)u(x, 0) = u_0(x) on an infinite or periodic domain, the exact analytic solution is simply the initial profile translated by a distance ctc t:

u(x,t)=u0(xct)u(x,t) = u_0(x - c t)

For a periodic domain, the coordinate xctx - c t is wrapped around the domain boundaries.

Discretization

To solve this numerically with the Finite Difference Method (FDM), we discretize the spatial domain into uniform cells of size Δx\Delta x, where xi=iΔxx_i = i \Delta x. Time is discretized into levels tn=nΔtt^n = n \Delta t.

We use a forward difference for the time derivative and a backward difference for the spatial derivative. The spatial difference is called an upwind scheme because we look "upwind" (in the direction the flow is coming from) to determine the gradient. Since c>0c > 0, information travels from left to right, so node ii depends on node i1i-1.

uin+1uinΔt+cuinui1nΔx=0\frac{u_i^{n+1} - u_i^n}{\Delta t} + c \frac{u_i^n - u_{i-1}^n}{\Delta x} = 0

Rearranging to explicitly solve for the new time level uin+1u_i^{n+1}:

uin+1=uinσ(uinui1n)u_i^{n+1} = u_i^n - \sigma (u_i^n - u_{i-1}^n)

Where σ\sigma is the Courant number.

Understanding the CFL number

The Courant-Friedrichs-Lewy (CFL) number defines the ratio of the physical wave speed to the numerical grid speed:

CFL=cΔtΔx\text{CFL} = \frac{c \Delta t}{\Delta x}

In this equation, cc is the convection speed, Δt\Delta t is the time step, and Δx\Delta x is the grid spacing. In simple terms, the CFL number measures how far information travels during one time step relative to the size of one grid cell.

Physically and numerically, if the CFL number is too large, the numerical method attempts to move information farther than the computational stencil can support in a single step. For the explicit first-order upwind scheme used here, the stable range is approximately 0CFL10 \le \text{CFL} \le 1 for positive cc.

Three Instructive Cases

  • CFL = 0.8: The scheme is stable, and the solution is transported correctly. However, the pulse becomes smeared because the first-order upwind method introduces numerical diffusion.
  • CFL = 1.0: The scheme is stable. For this specific periodic grid setup, the numerical pulse aligns very closely with the analytic shifted pulse. This is a special case where one time step corresponds exactly to one grid-cell shift. Do not overgeneralize this behavior to all schemes or all problems.
  • CFL = 1.2: This violates the CFL stability limit. As a result, the solution develops growing oscillations. In our visual results, the plotted values are clipped only to keep the diagram readable, while the printed error/min/max console values reveal the actual magnitude of the instability.

What the plot teaches

The CFL number controls both the stability and accuracy behavior of the numerical solution. Satisfying CFL1\text{CFL} \le 1 is strictly necessary for this scheme. Even within the stable range, a smaller stable CFL can still be more diffusive after many steps. The CFL = 1.2 case is intentionally included as a failure case to demonstrate what happens when the stability limit is exceeded.

Implementation Steps

  1. Setup parameters: Choose the domain size, grid resolution, wave speed, final time, and a stable CFL number.
  2. Initial Condition: Define the starting profile u0(x)u_0(x).
  3. Grid calculations: Compute Δx\Delta x and Δt\Delta t.
  4. Time Marching: Loop over time steps using the upwind update formula.
  5. Boundary Conditions: Apply periodic boundaries so the profile wraps around.
  6. Compare: Compute the analytic solution at the final time.
  7. Evaluate: Calculate the L2 error and plot the comparison.

Discussion of Results

The continuous solution translates without any shape change. The numerical upwind solution also translates, but it becomes smeared over time due to numerical diffusion. This is an artifact of the first-order truncation error in the upwind scheme.

The provided Python code now produces a comprehensive CFL sweep figure saved at output/linear_convection_cfl_sweep.png. This plot overlays the results for CFL = 0.8, 1.0, and 1.2. The analytic solution is drawn using black x markers so that the overlaid CFL = 1.0 curve remains visible. Furthermore, the unstable curve (CFL = 1.2) is clipped in the plot to maintain readability, although the underlying numerical computation diverges wildly.

Educational Python Example

You can find the complete, runnable Python code for this numerical experiment in The Numerics Lab repository. It includes the explicit upwind update and plotting to visualize numerical diffusion against the analytic solution.

View Example on GitHub (PR Branch)

Note: If the code PR is not merged yet, the example folder will appear in the main branch after the code PR is merged.