Durran and Blossey 2012 theta advection

Durran & Blossey 2012 describe a compressible Boussinesq model that has a Charney--Phillips staggering and a simplistic centred differencing of theta. Equation 41 describes theta advection:

db/dt + <<u>^z delta_x b>^x + <<w>^z delta_z b>^z = 0

where <a>^x denotes a two-point average of a in the x direction, and delta_x is a two-point derivative in the x direction. It can easily be shown that, for a constant wind, this reduces to centred differencing.

I've implemented this advection scheme in python and performed the same horizontal advection test. I'm using the same three-stage Runge-Kutta time-stepping as OpenFOAM.

Here's the errors at t=10000s:

I'm using this result as a baseline for a Charney--Phillips OpenFOAM advection solver. We should be able to get results that are comparable to, or better than, Durran and Blossey.