2.1 Fourier-von Neumann Analysis of Kinematic Wave Shallow Water Finite Element Schemes
The continuity equation of linearized shallow water system is given by \(\frac{\partial h}{\partial t}+u\frac{\partial h}{\partial x}+v\frac{\partial h}{\partial y}=0\), where \(u={\alpha }_{x}\beta {h}^{\beta -1}\), \(v={\alpha }_{y}\beta {h}^{\beta -1}\), \({\alpha }_{x}=\sqrt{{S}_{0x}}/{n}_{m}\), \({\alpha }_{y}=\sqrt{{S}_{0y}}/{n}_{m}\), \(\beta =5/3\), S0x and S0y are the element slopes in x and y directions, u and v are the x and y components of velocity, h is the flow depth and nm is the Manning’s roughness coefficient (Anmala and Mohtar, 2011). Assuming the finite element solution by \(\stackrel{-}{h}\), we obtain the error as \(\epsilon (x,y,t)=h(x,y,t)-\stackrel{-}{h}(x,y,t)\). Introducing the interpolation functions defined by \(\stackrel{-}{h}=\sum _{j=1}^{N}{H}_{j}\left(t\right){N}_{j}(x,y)\) and the residual error is minimized as follows: \(\underset{0}{\overset{1}{\int }}\underset{0}{\overset{1}{\int }}\left[{\stackrel{-}{h}}_{t}+u{\stackrel{-}{h}}_{x}+v{\stackrel{-}{h}}_{y}\right]{N}_{j}\left(x,y\right)dxdy=0\). Upon simplification, the finite element system can be written as \(A{H}^{\text{'}}+{B}_{x}H+{B}_{y}H+f=0\). The matrices A, Bx, By are given by
$$A={\int }_{{{\Omega }}_{e}}^{}{N}^{T}Nd{{\Omega }}_{e}=\frac{ab}{36}\left[\begin{array}{c}\begin{array}{cc}4& 2\\ 2& 4\end{array}\begin{array}{cc} 1& 2\\ 2& 1\end{array}\\ \begin{array}{cc}1& 2 \\ 2& 1 \end{array}\begin{array}{cc} 4& 2\\ 2& 4\end{array}\end{array}\right]$$
$${B}_{x}={\int }_{{{\Omega }}_{e}}^{}{N}^{T}\frac{\partial N}{\partial x}d{{\Omega }}_{e}=\frac{b}{12}\left[\begin{array}{c}\begin{array}{cc}-2& 2\\ -2& 2\end{array}\begin{array}{cc} 1& -1\\ 1& -1\end{array}\\ \begin{array}{cc}-1& 1 \\ -1& 1 \end{array}\begin{array}{cc} 2& -2\\ 2& -2\end{array}\end{array}\right]$$
and\(\)\({B}_{y}={\int }_{{{\Omega }}_{e}}^{}{N}^{T}\frac{\partial N}{\partial y}d{{\Omega }}_{e}=\frac{a}{12}\left[\begin{array}{c}\begin{array}{cc}-2& -1\\ -1& -2\end{array}\begin{array}{cc} 1& 2\\ 1& 1\end{array}\\ \begin{array}{cc}-1& -2 \\ -2& -1 \end{array}\begin{array}{cc} 2& 1\\ 2& 2\end{array}\end{array}\right]\)
where \(a=\varDelta x\) and \(b=\varDelta y\) for consistent method using the bilinear reactangular interpolation functions \({N}_{i}=\left(1-\frac{x}{a}\right)\left(1-\frac{y}{b}\right),{N}_{j}=\left(\frac{x}{a}\right)\left(1-\frac{y}{b}\right),{N}_{k}=\left(\frac{x}{a}\right)\left(\frac{y}{b}\right),{N}_{m}=\left(\frac{y}{b}\right)\left(1-\frac{x}{a}\right)\). The element matrix A for the lumped finite element scheme is given by
\(A=\frac{ab}{4}\left[\begin{array}{c}\begin{array}{cc}1& 0\\ 0& 1\end{array}\begin{array}{cc} 0& 0\\ 0& 0\end{array}\\ \begin{array}{cc}0& 0 \\ 0& 0 \end{array}\begin{array}{cc} 1& 0\\ 0& 1\end{array}\end{array}\right]\) . Following a similar analysis for error equations represented by \(\frac{\partial \epsilon }{\partial t}+u\frac{\partial \epsilon }{\partial x}+v\frac{\partial \epsilon }{\partial y}=0\), we obtain the error element and error gradient matrices. Approximating the error by Fourier series given by \({\int }_{{{\Omega }}_{e}}^{}\epsilon \left(x,y,t\right)d{{\Omega }}_{e}=\sum _{n=-\infty }^{n=+\infty }\sum _{m=-\infty }^{m=+\infty }{e}^{i\left({\omega }_{m}+{\omega }_{n}\right)t}{e}^{i({\sigma }_{m}x+{\sigma }_{n}y)}\), where \({\omega }_{m}\)and \({\omega }_{n}\) are temporal wave numbers and \({\sigma }_{m}\)and \({\sigma }_{n}\)are spatial wave numbers in x and y directions. This series can be approximated using individual error components such as \({\epsilon }_{ij}^{n}\left(x,y,t\right)={e}^{at}{e}^{i\left({\sigma }_{x}\varDelta x+{\sigma }_{y}\varDelta y\right)}\), and \({\epsilon }_{ij}^{n-1}\left(x,y,t\right)={e}^{-a\varDelta t}{\epsilon }_{ij}^{n}\left(x,y,t\right)\). The error equation can be solved for the amplification factor using the temporal wave number a. The specifics of stability analysis can be found from Press et al. (1992).
2.2 Theory: Lumped vs. Consistent Finite element schemes
Lumped mass matrix approach gives lower estimates of maximum eigen-values than distributed mass matrix or consistent mass matrix approach. Therefore, the lumped mass matrix would allow a larger time-step or critical time-step than the consistent mass matrix approach. In other words, the stability of a numerical scheme is assured with a larger time-step using mass lumping of capacitance matrix (Anmala and Mohtar, 2011). The issues of convergence and stability of lumped and distributed mass matrix approaches for linear and nonlinear problems are discussed by Wood (1990). The reduction of eigen-values for first and second-order scalar equations could be seen from distributed mass matrix modeling approach to lumped mass matrix modeling approach. In the same context, the effects of mass lumping for first and second-order problems, computation of eigen-values using Newmark and θ-method are detailed by Wood (1990).
The lower and upper bounds determined using the coefficient method for consistent and lumped element and node solutions are compared with integer multiples of eigen-solutions of consistent approximation. The closeness between the bounding solutions for all schemes with eigen-solutions of consistent approximation shows that the square finite element behaves essentially like a node at a larger scale. Node-node comparisons of amplification factors imply the accuracy of consistent finite element schemes and the alternate lumped finite element schemes. The solutions of amplification factors imply that the linearization of mass matrix coefficients in capacitance matrix and introduction of positive definiteness into the capacitance matrix do not alter the schemes significantly and only helps in obtaining better convergence rate.
2.3 Bounds on Amplification factors
The lower and upper bounds on amplification factors using the coefficient method are discussed in Anmala and Mohtar (2011). The lumped mass matrix assumption gives only one bounding estimate on Courant number using the coefficient method for all of the methods: explicit, semi-implicit, and implicit. The eigenvalue analysis of the system of finite element equations of shallow water for all the schemes (explicit, semi-implicit, implicit) and methods (consistent, lumped, upwind) could be obtained from that of Anmala and Mohtar (2015).