Finite-Differencing the 1D Wave Equation

March a vibrating string forward with a leapfrog scheme and meet the Courant condition that keeps it stable.

The equation

The wave equation u_tt = c^2 u_xx models a string under tension. Unlike diffusion it is second order in time, so we need two time levels of history to start.

Leapfrog scheme

Centered second differences in both time and space give u_new[i] = 2u[i] - u_old[i] + C^2 (u[i+1] - 2u[i] + u[i-1]), where the Courant number C = c dt/dx. The first step uses the initial velocity to bootstrap u_old.

Courant condition

Stability requires C <= 1: a wave must not cross more than one grid cell per time step. At C = 1 the scheme is exact for the pure advection pieces; above 1 it blows up. This CFL limit governs every explicit hyperbolic solver.

import numpy as np
nx=201; dx=1/(nx-1); c=1.0; C=0.9; dt=C*dx/c
x=np.linspace(0,1,nx)
u=np.exp(-300*(x-0.5)**2); u[0]=u[-1]=0
uold=u.copy()                       # zero initial velocity
for n in range(400):
    unew=np.zeros_like(u)
    unew[1:-1]=2*u[1:-1]-uold[1:-1]+C*C*(u[2:]-2*u[1:-1]+u[:-2])
    uold,u=u,unew
print('energy roughly conserved; peak',round(u.max(),3))

Behaviour

A Gaussian pulse splits into two half-amplitude pulses travelling left and right at speed c, reflecting off the fixed ends with a sign flip. Unlike the heat equation nothing decays - the scheme should conserve energy up to small dispersion, which is a useful sanity check on your implementation.