Calculating fluxes in local Lax-Friedrichs FVM for shallow water problem

finite-volume-methodfluid dynamicsnumerical methodspartial differential equations

so I have a more programming and some applied maths background, but absolutely none in FVM or Fluid dynamics. I have to implement a solver using the local Lax-Friedrichs method, and so far kinda-so-good. But I got stuck in calculating the fluxes (F and G) in the following Finite Volume solver equation:

$$U_{i,j}^{n+1} – U_{i,j}^{n} = \frac{\Delta t}{\Delta x}[F(u(x_{i+\frac{1}{2}},y,t))-F(u(x_{i-\frac{1}{2}},y,t))] + \frac{\Delta t}{\Delta y}[G(u(x,y_{i+\frac{1}{2}},t))-G((x,y_{i-\frac{1}{2}},t))]$$

The original differential Equation(which neglects a lot of forces) is as follows:

$$\frac{\partial}{\partial t}\begin{bmatrix}h \\ hu \\ hv \end{bmatrix}+ \frac{\partial}{\partial x}\begin{bmatrix}hu \\ hu^2 + 0.5gh^2 \\ huv \end{bmatrix}+ \frac{\partial}{\partial y}\begin{bmatrix}hv \\ huv \\ hv^2 + 0.5gh^2\end{bmatrix} = \begin{bmatrix} 0 \\ -\frac{\partial (ghb)}{\partial x}\\ -\frac{\partial (ghb)}{\partial y}\end{bmatrix}$$

Where:
b is the depth at a point
h is the distance from the bottom of the body of water till the lowest point of a "wave"
hu is one wavelength of a propagating wave\

I hope I explained it thoroughly, to sum up, I sort off understand what needs to be done, but I am having troubles on how to calculate the values of the fluxes to be used in the FV formula.

Best Answer

Rewrite the PDE as follows: $$\frac{\partial}{\partial t}U+ \frac{\partial}{\partial x}F(U)+ \frac{\partial}{\partial y}G(U) = 0$$ where $$ U = \begin{bmatrix}h \\ hu \\ hv \end{bmatrix},\quad F(U) = \begin{bmatrix}hu \\ hu^2 + 0.5gh^2 + ghb \\ huv \end{bmatrix}, \quad G(U) = \begin{bmatrix}hv \\ huv \\ hv^2 + 0.5gh^2 + ghb\end{bmatrix} . $$ The 2D local Lax-Friedrichs (LLF) method can be derived from these expressions, see this document. The scheme reads $$ U_{i,j}^{n+1} = U_{i,j}^n - \frac{\Delta t}{\Delta x}\left(F_{i+\frac12,j}^n - F_{i-\frac12,j}^n\right) - \frac{\Delta t}{\Delta y}\left(G_{i,j+\frac12}^n - G_{i,j-\frac12}^n\right) $$ where the LLF numerical fluxes are given in the linked document p. 33.