First-order implicit non-linear kinematic wave

Starting equations

Continuity with lateral inflow $q$: $\dfrac{\partial A}{\partial t}+\dfrac{\partial Q}{\partial x}=q$. The kinematic approximation $S_f=S_0$ with Manning's equation gives the unique relation

$$A=\alpha Q^{\beta},\qquad \alpha=\left(\frac{nP^{2/3}}{S_0^{1/2}}\right)^{0.6},\quad \beta=0.6$$

so that the non-linear kinematic wave equation is

$$\frac{\partial(\alpha Q^\beta)}{\partial t}+\frac{\partial Q}{\partial x}=q$$

Grid and differences

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

The equation is written at the point $(i+1,\ n+1)$, using

• a backward difference in time: $\dfrac{\partial A}{\partial t}\approx\dfrac{A_{i+1}^{n+1}-A_{i+1}^n}{\Delta t}$
• a backward difference in space, evaluated at the new time level $n+1$: $\dfrac{\partial Q}{\partial x}\approx\dfrac{Q_{i+1}^{n+1}-Q_i^{n+1}}{\Delta x}$
• $q\approx\bar q$ (average lateral inflow per unit length over the step).

Both differences have errors $O(\Delta t)$ and $O(\Delta x)$, so the scheme is first-order accurate. Because the spatial difference uses new-time values it is implicit and unconditionally stable.

Finite difference equation

$$\frac{\alpha\left(Q_{i+1}^{n+1}\right)^{\beta}-\alpha\left(Q_{i+1}^{n}\right)^{\beta}}{\Delta t}+\frac{Q_{i+1}^{n+1}-Q_i^{n+1}}{\Delta x}=\bar q$$

Multiplying by $\Delta t$ and putting the unknown $Q_{i+1}^{n+1}$ on the left:

$$\frac{\Delta t}{\Delta x}\,Q_{i+1}^{n+1}+\alpha\left(Q_{i+1}^{n+1}\right)^{\beta}=\frac{\Delta t}{\Delta x}\,Q_i^{n+1}+\alpha\left(Q_{i+1}^{n}\right)^{\beta}+\Delta t\,\bar q$$

The right side contains only known values ($Q_i^{n+1}$ from the upstream point at the new time, $Q_{i+1}^n$ from the previous time).

Solution (Newton-Raphson)

Let $x=Q_{i+1}^{n+1}$ and define

$$f(x)=\frac{\Delta t}{\Delta x}x+\alpha x^\beta-\left[\frac{\Delta t}{\Delta x}Q_i^{n+1}+\alpha\left(Q_{i+1}^n\right)^\beta+\Delta t\,\bar q\right]=0$$ $$x_{k+1}=x_k-\frac{f(x_k)}{f'(x_k)},\qquad f'(x)=\frac{\Delta t}{\Delta x}+\alpha\beta x^{\beta-1}$$

Start with $x_0=Q_{i+1}^n$ and iterate until the change is negligible. The calculation proceeds from the upstream boundary ($i=1$) down the channel for each time step.