Chapter 4 ended with a cost: about 17 points per wavelength to hold a drum's modes within 5 cents, and a factor of 2.8 in work for every halving of the error beyond that. The standard response is a more accurate stencil, and this chapter measures three versions of that response against each other. The results are less obvious than the orders of accuracy suggest, and one of the three is a trap.

Three schemes#

All three are explicit, advance two time levels, and use the rim condition $u = 0$.

Second order, the scheme of Chapter 4: the five-point Laplacian and leapfrog in time. Stable for $\lambda \le 1/\sqrt{2}$.

Fourth order in space only. Replace the three-point second difference in each direction by the five-point one,

which is fourth-order accurate, and keep leapfrog in time. The wider stencil has a larger largest eigenvalue, $32/(3h^2)$ against $8/h^2$, so the stable time step shrinks, to $\lambda \le \sqrt{3/8} \approx 0.612$. This is the scheme one arrives at by improving the part of the discretisation that is easiest to improve.

Fourth order in space and time. Leapfrog's leading time error is $(\Delta t^2/12)\,u_{tttt}$, and by the wave equation $u_{tttt} = c^4\nabla^4 u$, so it can be cancelled by adding $(\Delta t^4 c^4/12)\,\nabla^4 u$ to the update. Because that term already carries a factor $\Delta t^4$, a second-order approximation of $\nabla^4$ is enough, and the five-point Laplacian applied twice serves. The result is fourth order in both $h$ and $\Delta t$ for one extra stencil pass per step (Dablain 1986; it is often called the modified-equation or Lax–Wendroff approach). A short calculation with the symbols shows that its stability limit is $\lambda \le 1/\sqrt{2}$, the same as the second-order scheme's: the extra term gives back the time step that the wider stencil took away.

The wide stencils need values beyond the rim. With $u = 0$ on the rim, the wave equation forces $u_{tt} = 0$ there too, hence $\nabla^2 u = 0$, and repeating the argument shows that every even derivative vanishes on the rim, so the solution continues across it as an odd reflection. The script fills its ghost nodes that way. For the modes of the square this is exact, and the sine vectors remain exact eigenvectors of all three schemes, so pitch errors can be measured one mode at a time as before.

The measurement#

Three test modes, $(4,3)$, $(9,2)$ and $(14,1)$, cover the range from well resolved to barely resolved, and the last is within $4°$ of the $x$ axis, which Chapter 4 showed is the worst direction. Each is simulated alone for twelve of its periods, on grids from 16 to 96 intervals, with each scheme at a Courant number just inside its own stability limit: 0.7, 0.6 and 0.7. The error plotted is the worst of the three; the time is that of the run of the lowest mode, the longest of the three.

The second-order scheme behaves as it did in Chapter 4: its worst error falls from 409 cents at $N = 16$ to 7.7 at $N = 96$, with an observed order of 2.03 at the end.

The fully fourth-order scheme converges at order 4.0 over the whole ladder, from 245 cents to 0.19. At $N = 48$ its worst error is 3.0 cents, which the second-order scheme does not reach on any grid here. It is also more expensive per step, by a factor of about four and a half in this implementation, which is more than the extra stencil pass alone would suggest; part of that is the ghost-node handling that the wider stencils need, and a more careful implementation would narrow it. Even so, it wins at every accuracy that matters. Interpolating along the fourth-order curve, and extrapolating the second-order one beyond the grids measured at its observed order, 5 cents costs about 0.1 s with the fourth-order scheme and about 0.5 s with the second-order one, and 1 cent about 0.4 s against about 6 s. The advantage grows as the tolerance tightens, which is what the orders predict: halving the error costs a factor $2^{3/4} \approx 1.7$ in work at fourth order and $2^{3/2} \approx 2.8$ at second.

The scheme that is fourth order in space only is the trap. At $N = 24$ its worst error is 23 cents, much better than the second-order scheme's 146, and on a coarse grid it looks like a good idea. But leapfrog's time error is still there, now uncancelled by anything, and it has the opposite sign to the spatial error: the spatial error makes waves slow and the time error makes them fast. On the coarsest grids the two partly cancel, which produces the erratic observed orders in the table, and on fine grids the time error is left on its own, so the scheme converges at second order after all. By $N = 96$ its worst error, 5.1 cents, is barely better than the plain scheme's 7.7, at two and a half times the cost, because of the smaller time step and the wider stencil. The hover read-outs on the plots give the per-mode errors: at $N = 96$ the lowest test mode is 0.7 cents sharp under this scheme and 0.1 cents flat under the plain one.

This is a general lesson rather than a peculiarity of the wave equation. A method is as accurate as its least accurate part, and improving the spatial discretisation of a time-dependent problem without the temporal one improves only the constant in front of the error, and not always that.

Listening#

The same struck drum as in Chapter 4, on the same grid, with each scheme.

On a grid of $N = 32$ the relative rms error of the recording is 0.69 with the second-order scheme, 0.81 with fourth order in space only, and 0.37 with the fully fourth-order scheme; on $N = 24$ the fourth-order scheme's 0.46 is better than the second-order scheme's 0.59 on $N = 48$, for less computation. The waveform error is dominated, as before, by the highest modes the blow excites, which are poorly resolved on every grid here, so these numbers flatter none of the schemes. Listen for the upper partials of the ring rather than for its fundamental.

What this chapter establishes#

For the wave equation, where accuracy is lost to phase error accumulated over many periods, raising the order pays, and it pays more the tighter the tolerance. It pays only if the time discretisation is raised with the spatial one. The modified-equation scheme does this for the price of one more stencil pass, keeps the second-order scheme's time step, and here is roughly five times cheaper at 5 cents and fifteen times cheaper at 1 cent.

We emphasise that, on the square, every one of these schemes loses to the modal sum of Chapter 1, which is exact and costs one oscillator per mode. The square is here because it lets errors be measured exactly. Chapter 6 moves to a domain where the grid no longer fits the boundary, which is where the rest of the error comes from in practice.