Shock Capturing Methods for Shallow Water Equations

Question

I am looking for some help finding a numerical solution to the shallow water equations:

$\partial_tu+\partial_x(u^2/2+g\eta)=0$

$\partial_t \eta+\partial_x(u\eta)=0$.

where $u$ is the depth averaged horizontal velocity, $\eta$ is the depth, and $g$ is a constant.

Using mathematica's NDSolve I quickly develop a compressive shock (and consequentially the solution begins to ring and diverge) and I need a new solver to capture it accurately. It is my current understanding that I can capture the shock in at least three ways:

(a) I can capture the shock to arbitrary precision and guarantee uniqueness by using an artificial viscosity term. For this I have two questions:

1. What type of viscosity term should I use ($\partial_x^{2n} u$)? And what are the advantages and disadvantages.

2. Should I build a finite difference scheme and if so, what?

(b) I can capture the shock by having a flag that triggers when the gradient is large which then uses the Rankine–Hugoniot conditions to evolve the shock.

1. How do I decide on a trigger for the shock?

2. Should I build a finite difference scheme and if so, what?

(c) I can use a finite difference scheme that has enough numerical diffusion to not run into issues.

1. What scheme and how do I analyze the numerical diffusion contribution.

I don't know where to start. If you have any insight (maybe a solver on GitHub or a critique of that question itself) please share! Any advice on resources or methods would be helpful.

Answer 1

Educational resources

The shallow water equations are one of the simplest and most important systems of hyperbolic PDEs, and are frequently used in textbooks as an introductory example. There are several books that cover them from a theoretical viewpoint, but since you are primarily interested in numerical methods, I recommend LeVeque's finite volume book as well as Riemann Problems and Jupyter Solutions (full disclosure: I'm an author of the latter book). Toro's book is also excellent but out of print.

If your goal is to write your own solver (perhaps for educational purposes), then these books are a good place to learn what you'll need. If your goal is to solve some application problem, then you might prefer to use an existing code.

Solvers

The Clawpack software (disclosure: I'm one of the authors and principal maintainers) is often used for solving the shallow water equations and there are several example setups for that. It's probably easiest to get started using PyClaw. An example of a solution of the 2D shallow water equations solved in PyClaw is here, and more examples (both 1D and 2D, with or without varying bathymetry) can be found in the PyClaw examples directory (all the folders starting with "shallow").

For more advanced geophysical applications, you might want to use the GeoClaw code (also part of Clawpack), which is designed for modeling tsunamis, storm surges, etc.

There are many other solvers that can be found on the web, but I don't have direct experience with them.