Second-order scheme, linear kinematic wave

Linear kinematic wave equation

$$\frac{\partial Q}{\partial t}+c\,\frac{\partial Q}{\partial x}=c\,q$$

where $c=dQ/dA$ is taken as constant over the step and $q$ is the lateral inflow per unit length (zero if none). To get second-order accuracy, the derivatives are centred at the middle of the cell $(i+\tfrac12,\ n+\tfrac12)$ using all four corner points.

 t
 n+1 :  (i,n+1) o-------o (i+1,n+1)  <- unknown
                |   *   |
 n   :  (i,n)   o-------o (i+1,n)
              i         i+1       --> x

Central differences about the cell centre

$$\left.\frac{\partial Q}{\partial t}\right|_{i+\frac12}^{n+\frac12}\approx\frac{(Q_i^{n+1}-Q_i^n)+(Q_{i+1}^{n+1}-Q_{i+1}^n)}{2\Delta t},\qquad \left.\frac{\partial Q}{\partial x}\right|_{i+\frac12}^{n+\frac12}\approx\frac{(Q_{i+1}^{n+1}-Q_i^{n+1})+(Q_{i+1}^n-Q_i^n)}{2\Delta x}$$

Expanding each of the four values in a Taylor series about the centre, the first-derivative terms are retained and all second-derivative terms cancel, so the truncation errors are $O(\Delta t^2)$ and $O(\Delta x^2)$. The lateral inflow is also averaged:

$$q\approx\bar q=\tfrac14\left(q_i^n+q_{i+1}^n+q_i^{n+1}+q_{i+1}^{n+1}\right)$$

Finite difference equation

$$\frac{Q_i^{n+1}+Q_{i+1}^{n+1}-Q_i^n-Q_{i+1}^n}{2\Delta t}+c\,\frac{Q_{i+1}^{n+1}+Q_{i+1}^n-Q_i^{n+1}-Q_i^n}{2\Delta x}=c\,\bar q$$

Multiply by $2\Delta t$ and let $C=c\Delta t/\Delta x$:

$$(1+C)\,Q_{i+1}^{n+1}=(1+C)\,Q_i^n+(1-C)\left(Q_{i+1}^n-Q_i^{n+1}\right)+2c\,\Delta t\,\bar q$$

Result: discharge at the unknown time and location

$$\boxed{\,Q_{i+1}^{n+1}=Q_i^n+\frac{1-C}{1+C}\left(Q_{i+1}^n-Q_i^{n+1}\right)+\frac{2C\,\Delta x}{1+C}\,\bar q\,}$$

where the last term is for lateral inflow. Check of the coefficient: $2c\Delta t\bar q/(1+C)=2C\Delta x\bar q/(1+C)$.

Remarks

• $Q_i^n$, $Q_{i+1}^n$ are known from the previous time level and $Q_i^{n+1}$ from the upstream boundary (or the previous point), so $Q_{i+1}^{n+1}$ follows directly: march downstream along $i$ for each time step.
• Stability: the scheme is unconditionally stable. Accuracy is best for $C\approx1$; for $C\gg1$ or $C\ll1$ it can show numerical oscillations (dispersion).
• If $C=1$, the formula gives $Q_{i+1}^{n+1}=Q_i^n$ (exact translation of the wave).