Support other models

[JoshuaLampert]

It would be nice to support other dispersive wave equations. These include especially the equations in H. Ranocha, D. Mitsotakis, D. Ketcheson, A Broad Class of Conservative Numerical Methods for Dispersive Wave Equations (2021) and the corresponding discretizations using SBP operators developed there. These include the following scalar equations with constant bottom:

• BBM (Add BBMEquation1D #150)
• Fornberg-Whitham
• Camassa-Holm
• Degasperis-Procesi
• Holm-Hone

Further models are discussed in J. Giesselmann, H. Ranocha, Convergence of hyperbolic approximations to higher-order PDEs for smooth solutions (2025):

• Korteweg-de Vries-Burgers
• Gardner
• Kawahara
• Generalized Kawahara
• Linear bi-harmonic
• Kuramoto-Sivashinsky

and their hyperbolizations.
There are also the

• Serre-Green-Naghdi equations in different formulations (H. Ranocha, M. Ricchiuto, Structure-preserving approximations of the Serre-Green-Naghdi equations in standard and hyperbolic form (2024)), see Roadmap for Serre-Green-Naghdi #129.

Of course the classical

• Korteweg-de Vries (KdV) equation (there is some code at https://github.com/abhibsws/2024_kdvh_RR) (Implementation of the 1D KdV equation.  #198, Roadmap KdV eq. #202 )

would be of interest, too as well as the

• hyperbolic approximation of the Korteweg-de Vries (KdVH) equation (A. Biswas, D. Ketcheson, H. Ranocha, J. Schütz, Traveling-wave solutions and structure-preserving numerical methods for a hyperbolic approximation of the Korteweg-de Vries equation (2024)). There is some code at https://github.com/abhibsws/2024_kdvh_RR
• hyperbolic approximation of the BBM (BBMH) equation (S. Bleecke, A. Biswas, D. Ketcheson,
  H. Ranocha, J. Schütz, Asymptotic-preserving and energy-conserving methods for a hyperbolic approximation of the BBM equation (2025)). There is code at https://github.com/sbleecke/2025_bbmh.

There is also

• a hyperbolic approximation of the Camassa-Holm equation (CHH) is given in A. Biswas, D. Ketcheson, Approximation of arbitrarily high-order PDEs by first-order hyperbolic relaxation (2025)
• the nonlinear Schrödinger (NLS) equation and its hyperbolizaton from H. Ranocha, D. Ketcheson, High-order mass- and energy-conserving methods for the nonlinear Schrödinger equation and its hyperbolization (2025) with code at https://github.com/ranocha/2025_nls

Implementing these semidiscretizations should be pretty straightforward within DispersiveShallowWater.jl.
Another interesting model to look at (requiring to develop energy preserving semidiscretizations using SBP operators first) would be the

• equations by Escalante and Morales de Luna (Escalante and Morales de Luna, A General Non-hydrostatic Hyperbolic Formulation for Boussinesq Dispersive Shallow Flows and Its Numerical Approximation (2020))

[JoshuaLampert]

I'm unsure what the best approach would be to add the KdV equation. So far we have implemented the equations in dimensional form, e.g. the BBM equation as

$$ \eta_t + \sqrt{gD}\eta_x + 3/2\sqrt{g/D}\eta\eta_x - 1/6D^2\eta_{xxt} = 0. $$

In the literature you often find the more abstracted version

$$ u_t + u_x + uu_x - u_{xxt} = 0. $$

That's fine because these two versions are the same except for constants, which means you can (basically) solve the abstracted form, e.g., by picking $g = D = 1$.
However, for the KdV equation this is a bit different. I would like to support the dimensional form of the equation, such that $\eta$ has the same (physical) meaning across equations in order to be able to compare the solutions of different equations with each other. This means the KdV equation would look like (see also #150 (comment)):

$$ \eta_t + \sqrt{gD}\eta_x + 3/2\sqrt{g/D}\eta\eta_x + 1/6\sqrt{gD}D^2\eta_{xxx} = 0. $$

The form usually dealt with in the literature looks like (sometimes with an additional factor -6 in front of the $uu_x$ term)

$$ u_t + uu_x + u_{xxx} = 0, $$

which means it does not contain a term including only one spatial derivative of the solution. This means if we implement the dimensional form, it is not so straightforward to solve the common non-dimensional version usually used in papers, which is not very nice.
Do you see a way out of this dilemma, @ranocha? Could we maybe provide a nice interface computing $\eta$ from $u$ and vice versa, which allows solving both versions? I haven't thought this through though.

[ranocha]

What about the transformation $u = \eta - 1$ or something like that?

[JoshuaLampert]

So with $u = \eta/D + 2/3$, we would get

$$ u_t + 3/2\sqrt{gD}uu_x + 1/6\sqrt{gD}D^2u_{xxx} = 0, $$

which would end up being $u_t + uu_x + u_{xxx} = 0$ for $D = 3$, $g = 4/27$.
Could we use prim2phys for the transformation, i.e. the primitive variable is $u$ in this case, which is also used within the discretization and we compute $\eta$ from $u$ by calling eta = prim2phys(q, equations)? Then we would also need to specialize waterheight_total and so on.
Maybe this could work out, but I would need to test it.

[JoshuaLampert]

As initial condition one could either set directly u of course or set the physical waterheight and then transforming it with a function phys2prim.

[ranocha]

Sounds good to me. To clarify - do you want to use $u$ or $\eta$ as a variable in the solver? Both options are fine with me - but remember that we used $\eta$ so far.

[JoshuaLampert]

I wanted to propose to use $u$ in the solver. That's why we would need to define specializations for waterheight_total etc. The main reason for why I propose to use $u$ as variable is that most discretizations are available for this formulation and therefore we would not need to translate the discretizations first. In other words, we always use the primitive variables in the solver.

[ranocha]

👍 We should clarify this in the docstrings but otherwise it's fine with me.

[JoshuaLampert]

Yes, definitely 👍

[JoshuaLampert]

After a discussion with @cwittens, we agreed that it's better to implement the dimensional variant of the KdV equation to be consistent with the other equations. This requires to translate the semidiscretization and quadratic invariant into the dimensional form, but this should be doable. To be able to compare it with the non-dimensional form, we could add a function prim2nondim (and the inverse function nondim2prim), which computes $$u$$ from $$\eta$$ by $$u = \eta/D + 2/3$$. With this transformation one should get the same results either running the implemented dimensional form and transforming the result or directly running the non-dimensional form. We should check whether that's indeed the case (e.g., by comparing to https://github.com/abhibsws/2024_kdvh_RR) once this is implemented.

[ranocha]

I agree 👍