Two-layer Shallow Water Equations Research Articles

An Entropy Stable Discontinuous Galerkin Method for the Two-Layer Shallow Water Equations on Curvilinear Meshes

We present an entropy stable nodal discontinuous Galerkin spectral element method (DGSEM) for the two-layer shallow water equations on two dimensional curvilinear meshes. We mimic the continuous entropy analysis on the semi-discrete level with the DGSEM constructed on Legendre–Gauss–Lobatto (LGL) nodes. The use of LGL nodes endows the collocated nodal DGSEM with the summation-by-parts property that is key in the discrete analysis. The approximation exploits an equivalent flux differencing formulation for the volume contributions, which generate an entropy conservative split-form of the governing equations. A specific combination of a numerical surface flux and discretization of the nonconservative terms is then applied to obtain a high-order path-conservative scheme that is entropy conservative. Furthermore, we find that this combination yields an analogous discretization for the pressure and nonconservative terms such that the numerical method is well-balanced for discontinuous bathymetry on curvilinear domains. Dissipation is added at the interfaces to create an entropy stable approximation that satisfies the second law of thermodynamics in the discrete case, while maintaining the well-balanced property. We conclude with verification of the theoretical findings through numerical tests and demonstrate results about convergence, entropy stability and well-balancedness of the scheme.

A robust numerical scheme based on auxiliary interface variables and monotone-preserving reconstructions for two-layer shallow water equations with wet–dry fronts

This work aims to propose a well-balanced and positivity-preserving numerical scheme for the two-dimensional two-layer shallow water equations with a nonflat bottom topography on uniform meshes. One difficulty is to preserve the water at rest when the computational domain contains wet–dry fronts. It is highly nontrivial. Another difficulty in solving the system is to deal with the nonconservative product term which cannot be defined in the distributional sense. To overcome these issues, we define auxiliary interface variables motivated by the generalized hydrostatic reconstruction (Castro et al., 2007) for computing Riemann states to obtain consistent numerical fluxes across the cell interface. When the water is at rest, the numerical flux should be continuous across the cell interface and the nonconservative product term is equivalent to the gradient of physical fluxes. In order to maintain steady states at a discrete level, we discretize the integral of the source term in the upwind sense. One of the main contributions of this work is to introduce a monotone-preserving reconstruction for constructing non-oscillatory piecewise linear polynomials to improve the accuracy in space. The current scheme can guarantee both layer heights to be nonnegative and maintain the water at rest even if the computational domain contains wet–dry fronts. Several numerical results indicate that the current scheme is stable in solving the two-dimensional two-layer shallow water equations with wet–dry fronts. In particular, the current scheme is stable when computing complicated states which contain many shocks and will evolve to a turbulent state.

Two-Layer Shallow Water Equations with Momentum Conservative Scheme for Wave Propagation Simulation

In this paper, we discuss the implementation of momentum conservative scheme to shallow water equations (SWE). In shallow water model, the hydrodynamic pressure of the water is neglected. Here, the numerical calculation of mass and momentum conservation was applied on a staggered grid domain. The vertical interval was divided into two parts which made the computation quite efficient and accurate. Our focus is on the performance of the numerical scheme in simulating wave propagation and run-up phenomena, where the main challenge is to calculate the wave speed accurately and to count the non-linear term of the model. Here we also considered the wet and dry conditions of the topography. Three benchmark tests were picked out to validate the numerical scheme. A simulation of standing wave was carried out; the results were compared to the linear analytical solution and show a good fit. In addition, a simulation of harmonic wave propagation on a sloping beach was conducted, and the results closely align with the expected values from exact solution. Finally, we carried out a simulation of solitary wave with a sloping topography; and the results were compared to laboratory data. A good agreement was observed between the simulation results and experimental measurements.

A high-order domain preserving DG method for the two-layer shallow water equations

The bilayer shallow water wave equations in one-dimensional space are considered in this paper. The equations admit two groups of characteristic velocities, which are the first-order approximation of the eigenvalues. Due to the numerical instability, the characteristic velocities may become complex, and thus the system is not hyperbolic and yields to the so-called Kelvin–Helmholtz instability at the interface of the two layers. To overcome this issue, an invariant domain preserving DG method is presented for the bilayer shallow water wave equations. The proposed method is high-order accurate, conservative and can keep the characteristic velocities being real provided that the initial characteristic velocities are real. Therefore, the Kelvin–Helmholtz instability at the interface can be avoided. Representative numerical examples are chose to demonstrate the performance of the proposed method.

High-order compact gas-kinetic scheme for two-layer shallow water equations on unstructured mesh

For the two-layer shallow water equations, a high-order compact gas-kinetic scheme (GKS) on triangular mesh is proposed. The two-layer shallow water equations have complex source terms in comparison with the single layer equations. The main focus of this study is to construct a time-accurate evolution solution at a cell interface and to design a well-balanced scheme. The evolution model at a cell interface provides not only the numerical fluxes, but also the flow variables. In the development of the well-balanced scheme, the time-dependent flow variables at the cell interfaces surrounding the cell can be utilized to update the cell-averaged gradients. These time-evolving gradients are subsequently applied in the discretization of source terms inside each control volume. Based on the cell-averaged flow variables and their gradients, high-order initial data reconstruction can be achieved with compact stencils. The compact high-order GKS has advantages to simulate the flow evolution in complex domain covered by unstructured mesh. Many test cases are used to validate the accuracy and robustness of the scheme for the two-layer shallow water equations.

Preserving stationary discontinuities in two-layer shallow water equations with a novel well-balanced approach

Abstract This paper proposes a novel energy-balanced numerical scheme for the two-layer shallow water equations (2LSWEs) that accurately captures internal hydraulic jumps without introducing spurious oscillations. The proposed scheme overcomes the problem of post-shock oscillations in the 2LSWE, a phenomenon commonly observed in numerical solutions of non-linear hyperbolic systems when shock-capturing schemes are used. The approach involves reconstructing the internal momentum equation of 2LSWEs using the correct Hugoniot curve via a set of shock wave fixes originally developed for single-layer shallow water equations. The scheme successfully preserves all stationary solutions, making it highly suitable for simulations of real-life scenarios involving small perturbations of these conditions.

A hybrid 1D-2D Lagrangian solver with moving coupling to simulate dam-break flow

In dam-break problems, two-dimensional (2D) models on the vertical plane appropriately capture the non-hydrostatic pressure distribution but with a high computational cost. This study presents a hybrid Lagrangian solver of one-dimensional (1D) shallow water flow and a 2D large-eddy simulation with a sub-particle scale to simulate dam-break flows over frictional beds. The proposed 1D-2D hybrid scheme is a useful trade-off between the accuracy of the 2D models on the vertical plane and the computational efficiency of 1D shallow water equations models. The configuration of sub-domains is established so that the 1D solves the regions with dominantly hydrostatic pressure, while the 2D solver accounts for the non-hydrostaticity. A moving coupling method is introduced based on two-layer shallow water equations with a uniform velocity distribution that constructs two sub-domains, including inner and outer zones in the section coupling model (SCM), or upper and lower layers within the layer coupling model (LCM). Then, the region with high vertical accelerations is substituted by the 2D model in the SCM or mixed SCM-LCM in which the overlapping zone and boundary interface move with time, based on the 1D flow feature. Also, a new method is proposed to deal with irregular boundary interfaces in the LCM while satisfying the momentum transfer. The present hybrid model is advantageous as it meets mass conservation, and there is no need for typical inflow/outflow boundary conditions at the coupling interfaces for simulating dam-break waves. Moving particle simulation (MPS) is utilized to solve the Lagrangian equations, which is advantageous in dealing with free-surface problems with large deformation. The model validation reveals that the hybrid results of water depths are in very good agreement with observed data and those obtained from the 2D model for selected dam-break flows with the ratio of initial depths downstream to upstream up to 0.4. In this respect, the average root-mean-square-error (RMSE) is computed at 8.09 mm in the SCM, 8.16 mm in the mixed SCM-LCM, and 7.48 mm in the pure 2D model, while both hybrid models save the computational time by 55% compared to the 2D model. Moreover, the mixed SCM-LCM preserves the nonlinearity of the wavefront over long times, whereas the SCM converts the wavefront into a shockwave considering the same number of fluid particles.

On well-balanced implicit-explicit Lagrange-projection schemes for two-layer shallow water equations

This works concerns the study of well-balanced Lagrange-projection schemes applied to the two-layer shallow water system. In particular, a formulation of the mathematical model in Lagrangian coordinates is proposed. The HLL method is then applied to a simplified version of the resulting Lagrangian system. Furthermore, based on the acoustic-transport splitting interpretation, another approximate Riemann solver for the acoustic-Lagrangian step is described. Both an explicit and an implicit-explicit method are proposed, where the latter can allow very fast simulations in sub-critical regimes. Finally, we show numerical simulations in which the outputs are compared with the IFCP method’s results.

Well-balanced path-conservative central-upwind schemes based on flux globalization

In this paper, we introduce a new approach for constructing robust well-balanced (WB) finite-volume methods for nonconservative one-dimensional hyperbolic systems of nonlinear partial differential equations. The WB property, namely, the ability of the scheme to exactly preserve physically relevant steady-state solutions is enforced using a flux globalization approach according to which a studied system is rewritten in an equivalent quasi-conservative form with global fluxes. To this end, one needs to incorporate nonconservative product terms into the global fluxes. The resulting system can then be solved using a Riemann-problem-solver-free central-upwind (CU) scheme. However, a straightforward integration of the nonconservative terms would result in a scheme capable of exactly preserving very simple smooth steady states only and failing to preserve discontinuous steady states naturally arising in the nonconservative models.In order to ameliorate the flux globalization based CU scheme, we evaluate the integrals of the nonconservative product terms using the technique introduced in [5], where a path-conservative central-upwind scheme (PCCU) was introduced. This results in a new flux globalization based WB PCCU scheme, which is much more accurate and robust than both the original PCCU scheme and the straightforward flux globalization based CU scheme. This is illustrated using two nonconservative systems: a system describing fluid flows in nozzles with variable cross-sections and the two-layer shallow water equations. We demonstrate superiority of the proposed flux globalization based WB PCCU scheme on a number of challenging examples.

Fifth-order A-WENO schemes based on the path-conservative central-upwind method

We develop fifth-order A-WENO finite-difference schemes based on the path-conservative central-upwind method for nonconservative one- and two-dimensional hyperbolic systems of nonlinear PDEs. The main challenges in development of accurate and robust numerical methods for the studied systems come from the presence of nonconservative products. Semi-discrete second-order finite-volume path-conservative central-upwind (PCCU) schemes recently proposed in Castro Díaz et al. (2019) [8] provide one with a reliable Riemann-problem-solver-free numerical method for nonconservative hyperbolic system. In this paper, we extend the PCCU schemes to the fifth-order of accuracy in the framework of A-WENO finite-difference schemes.We apply the developed schemes to the two-layer shallow water equations. We ensure that the developed schemes are well-balanced in the sense that they are capable of exactly preserving “lake-at-rest” steady states. We illustrate the performance of the new fifth-order schemes on a number of one- and two-dimensional examples, where one can clearly see that the proposed fifth-order schemes clearly outperform their second-order counterparts.

A Nonintrusive Reduced-Order Model for Uncertainty Quantification in Numerical Solution of One-Dimensional Free-Surface Water Flows Over Stochastic Beds

Free-surface water flows over stochastic beds are complex due to the uncertainties in topography profiles being highly heterogeneous and imprecisely measured. In this study, the propagation and influence of several uncertainty parameters are quantified in a class of numerical methods for one-dimensional free-surface flows. The governing equations consist of both single-layer and two-layer shallow water equations on either flat or nonflat topography. For this purpose, the free-surface profiles are computed for different realizations of the random variables when the bed is excited with sources whose statistics are well defined. Many research studies have been dedicated to the development of numerical methods to achieve some order of accuracy in free-surface flows. However, little concern was given to examine the performance of these numerical methods in the presence of uncertainty. This work addresses this specific area in computational hydraulics with regards to the uncertainty generated from bathymetric forces. As numerical solvers for the one-dimensional shallow water equations, we implement four finite volume methods. To reduce the required number of samples for uncertainty quantification, we combine the proper orthogonal decomposition method with the polynomial chaos expansions for efficient uncertainty quantifications of complex hydraulic problems with large number of random variables. Numerical results are shown for several test examples including dam-break problems for single-layer and two-layer shallow water flows. The problem of flow exchange through the Strait of Gibraltar is also solved in this study. The obtained results demonstrate that in some hydraulic applications, a highly accurate numerical method yields an increase in its uncertainty and makes it very demanding to use in an operational manner with measured data from the field. On the other hand, when the complexity of physics increases, these highly accurate numerical methods display less uncertainty compared to the low accurate methods.

Well-balanced numerical resolution of the two-layer shallow water equations under rigid-lid with wet–dry fronts

In this paper, we present a well-balanced numerical scheme for the frictional two-layer shallow water equations (2LSWE) over variable bottom topography and under a rigid-lid (RL) to simulate internal waves propagating over wet and dry areas. Following the idea of Liang and Borthwick (2009) for the derivation of a pre-balanced formulation for one-layer flows, we derive a new formulation of the 2LSWE-RL using the interface elevation above datum and horizontal momentum as conservative variables. This new formulation mathematically balances the flux gradient and source terms so that the lake-at-rest steady state is automatically preserved in wet-bed applications. A proper discretization of the slope source term is adopted to produce well-balanced solutions in dry-bed areas where the lower layer depth vanishes. Using the Harten, Lax, and van Leer (HLL) approximate Riemann solver, we adopt a finite volume MUSCL-RK2 method to construct a second-order well-balanced numerical scheme ensuring the non-negativity of the layers depths with a special treatment of source terms. Finally, various numerical experiments are tested in both one and two-layer cases validating the properties of the numerical scheme.

A new well-balanced finite-volume scheme on unstructured triangular grids for two-dimensional two-layer shallow water flows with wet-dry fronts

In this study, a novel two-dimensional finite-volume method on unstructured triangular meshes for two-layer shallow water flows is developed. First, the one-dimensional relaxation approach developed in Abgrall and Karni (2009) [1] is extended to the currently studied two-dimensional two-layer shallow water equations to build an unconditionally hyperbolic system. Next, a two-dimensional finite-volume HLL scheme on triangular grids for such a non-conservative hyperbolic system is constructed. Special attention is paid to guaranteeing the well-balanced property even in the presence of the wet-dry fronts. To this end, techniques such as spatial reconstruction for the secondary variables, special cell classifications for the two-layer fluid system, and special local reconstruction of auxiliary surface variables are proposed. These special techniques are further combined with new discrete formulas for the nonconservative products and the bed-slope terms, and it thus leads to an exactly well-balanced numerical scheme for the studied two-layer system even in the presence of wet-dry fronts. The developed numerical model is second order accurate and preserves the nonnegative layer depths. The performances of the proposed numerical model are investigated by several numerical experiments.

Modeling Heavy-Gas Dispersion in Air with Two-Layer Shallow Water Equations

Computation of gas dispersal in urban places or hilly grounds requires a large amount of computer time when addressed with conventional multidimensional models. Those are usually based on two-phase flow or Navier-Stokes equations. Different classes of simplified models exist. Among them, two-layer shallow water models are interesting to address large-scale dispersion. Indeed, compared to conventional multidimensional approaches, 2D simulations are carried out to mimic 3D effects. The computational gain in CPU time is consequently expected to be tremendous. However, such models involve at least three fundamental difficulties. The first one is related to the lack of hyperbolicity of most existing formulations, yielding serious consequences regarding wave propagation. The second is related to the non-conservative terms in the momentum equations. Those terms account for interactions between fluid layers. Recently, these two difficulties have been addressed in Chiapolino and Saurel (2018) and an unconditional hyperbolic model has been proposed along with a Harten-Lax-van Leer (HLL) type Riemann solver dealing with the non-conservative terms. In the same reference, numerical experiments showed robustness and accuracy of the formulation. In the present paper, a third difficulty is addressed. It consists of the determination of appropriate drag effect formulation. Such effects also account for interactions between fluid layers and become of particular importance when dealing with heavy-gas dispersion. With this aim, the model is compared to laboratory experiments in the context of heavy gas dispersal in quiescent air. It is shown that the model accurately reproduces experimental results thanks to an appropriate drag force correlation. This function expresses drag effects between the heavy and light gas layers. It is determined thanks to various experimental configurations of dam-break test problems.

On the efficient implementation of PVM methods and simple Riemann solvers. Application to the Roe method for large hyperbolic systems

Polynomial Viscosity Matrix (PVM) methods can be considered as approximations of the Roe method in which the absolute value of the Roe matrix appearing in the numerical viscosity is replaced by the evaluation of the Roe matrix at a chosen polynomial that approximates the absolute value function. They are in principle cheaper than the Roe method since the computation and the inversion of the eigenvector matrix is not necessary. In this article, an efficient implementation of the PVM based on polynomials that interpolate the absolute value function at some points is presented. This implementation is based on the Newton form of the polynomials. Moreover, many numerical methods based on simple Riemann solvers (SRS) may be interpreted as PVM methods and thus this implementation can be also applied to them: the close relation between PVM methods and simple Riemann solvers is revisited here and new shorter proofs based on the classical interpolation theory are given. In particular, Roe method can be interpreted both as a SRS and as a PVM method so that the new implementation can be used. This implementation, that avoids the computation and the inversion of the eigenvector matrix, is called Newton Roe method. Newton Roe method yields the same numerical results of the standard Roe method, with less runtime for large PDE systems. Numerical results for two-layer Shallow Water Equations and Quadrature-Based Moment Equations show a significant speedup if the number of equations is large enough.

A High Order Central DG method of the Two-Layer Shallow Water Equations

A High Order Central DG method of the Two-Layer Shallow Water Equations

Re-evaluating efficiency of first-order numerical schemes for two-layer shallow water systems by considering different eigenvalue solutions

The efficiency of several first-order numerical schemes for two-layer shallow water equations (SWE) are evaluated here by considering different eigenvalue solutions. This study is a continuation of our previous work (Krvavica et al., 2018) in which we have proposed an efficient implementation of a Roe solver for two-layer SWE based on analytical expressions for eigenvalues and eigenvectors. In this work, the accuracy and computational cost of numerical, analytical, and approximated eigenvalue solvers are compared when implemented in Roe, Intermediate Field CaPturing (IFCP) and Polynomial Viscosity Matrix (PVM) schemes. Several numerical tests are performed to examine the overall efficiency of numerical schemes with different eigenvalue solvers when computing two-layer shallow-water flows. The results confirm that analytical eigenvalue solutions are much faster than numerical solvers, with a computational cost closer to approximate expressions. Consequently, the Roe scheme with analytical solutions to the eigenstructure is equally efficient as the IFCP scheme. On the other hand, IFCP and PVM schemes with analytical solutions to eigenvalues are found to be equally efficient as those with approximated expressions. Analytical eigenvalues show slightly better results when dealing with larger density differences between the layers.

A consistent reduction of the two-layer shallow-water equations to an accurate one-layer spreading model

The gravity-driven spreading of one fluid in contact with another fluid is of key importance to a range of topics. These phenomena are commonly described by the two-layer shallow-water equations (SWE). When one layer is significantly deeper than the other, it is common to approximate the system with the much simpler one-layer SWE. It has been assumed that this approximation is invalid near shocks, and one has applied additional front conditions to correct the shock speed. In this paper, we prove mathematically that an effective one-layer model can be derived from the two-layer equations that correctly capture the behavior of shocks and contact discontinuities without additional closure relations. The result shows that simplification to an effective one-layer model is justified mathematically and can be made without additional knowledge of the shock behavior. The shock speed in the proposed model is consistent with empirical models and identical to front conditions that have been found theoretically by von Kármán and Benjamin. This suggests that the breakdown of the SWE in the vicinity of shocks is less severe than previously thought. We further investigate the applicability of the SW framework to shocks by studying one-dimensional lock-exchange/-release. We derive expressions for the Froude number that are in good agreement with the widely employed expression by Benjamin. The equations are solved numerically to illustrate how quickly the proposed model converges to solutions of the full two-layer SWE. We also compare numerical results from the model with results from experiments and find good agreement.

A Two-Layer Hydrostatic-Reconstruction Method for High-Resolution Solving of the Two-Layer Shallow-Water Equations over Uneven Bed Topography

In this study, attention is focused on numerically solving the two-dimensional, two-layer nonlinear shallow-water equations (2LSWEs) over uneven bed topography. A two-layer hydrostatic-reconstruction method (2LHRM) is proposed for face value reconstructions. A numerical model is then developed based on the 2LHRM in the framework of finite-volume methods using the slope-limited centered (SLIC) approximate Riemann solver for flux evaluations. The validations against various benchmark tests show that the 2LHRM by working with the SLIC scheme predicts robust solutions around large gradients and is able to simulate lake-at-rest solutions without using any ad-hoc techniques.

The Influence of Topography on the Stability of the Norwegian Atlantic Current off Northern Norway

AbstractThe steep continental slope off the Lofoten–Vesterålen islands of northern Norway appears to be the source of the most intense mesoscale eddy field in all of the Nordic Seas. Here we use linearized two-layer shallow-water equations to study the stability of the Norwegian Atlantic Current in this region. The study extends previous works that used linearized quasigeostrophic vertical mode equations to examine the effects of bottom topography on baroclinic instability here. We find evidence of baroclinic instability in the stacked shallow-water model but also of barotropic instability that is associated with the upper-layer lateral shear. The calculations give an indication that growth rates of barotropic instability may be comparable to or larger than those of baroclinic instability over the steepest parts of the continental slope.