A Well-Balanced Positivity-Preserving Quasi-Lagrange Moving Mesh DG Method for the Shallow Water Equations

A high-order, well-balanced, positivity-preserving quasi-Lagrange moving mesh DG method is presented for the shallow water equations with non-flat bottom topography. The well-balance property is crucial to the ability of a scheme to simulate perturbation waves over the lake-at-rest steady state such as waves on a lake or tsunami waves in the deep ocean. The method combines a quasi-Lagrange moving mesh DG method, a hydrostatic reconstruction technique, and a change of unknown variables. The strategies in the use of slope limiting, positivity-preservation limiting, and change of variables to ensure the well-balance and positivity-preserving properties are discussed. Compared to rezoning-type methods, the current method treats mesh movement continuously in time and has the advantages that it does not need to interpolate flow variables from the old mesh to the new one and places no constraint for the choice of a update scheme for the bottom topography on the new mesh. A selection of one- and two-dimensional examples are presented to demonstrate the well-balance property, positivity preservation, and high-order accuracy of the method and its ability to adapt the mesh according to features in the flow and bottom topography.     

High-Order Positivity-Preserving Well-Balanced Discontinuous Galerkin Methods for Euler Equations with Gravitation on Unstructured Meshes

In this paper, we propose a high-order accurate discontinuous Galerkin (DG) method for the compressible Euler equations under gravitational fields on unstructured meshes. The scheme preserves a general hydrostatic equilibrium state and provably guarantees the positivity of density and pressure at the same time. Comparing with the work on the well-balanced scheme for Euler equations with gravitation on rectangular meshes, the extension to triangular meshes is conceptually plausible but highly nontrivial. We first introduce a special way to recover the equilibrium state and then design a group of novel variables at the interface of two adjacent cells, which plays an important role in the well-balanced and positivity-preserving properties. One main challenge is that the well-balanced schemes may not have the weak positivity property. In order to achieve the well-balanced and positivity-preserving properties simultaneously while maintaining high-order accuracy, we carefully design DG spatial discretization with well-balanced numerical fluxes and suitable source term approximation. For the ideal gas, we prove that the resulting well-balanced scheme, coupled with strong stability preserving time discretizations, satisfies a weak positivity property. A simple existing limiter can be applied to enforce the positivity-preserving property, without losing high-order accuracy and conservation. Extensive one- and two-dimensional numerical examples demonstrate the desired properties of the proposed scheme, as well as its high resolution and robustness.     

A Well-Balanced Gas Kinetic Scheme for Navier-Stokes Equations with Gravitational Potential

The hydrostatic equilibrium state is the consequence of the exact balance between hydrostatic pressure and external force. Standard finite volume cannot keep this balance exactly due to their unbalanced truncation errors. In this study, we introduce an auxiliary variable which becomes constant at isothermal hydrostatic equilibria and propose a well-balanced gas kinetic scheme for the Navier-Stokes equations. Through reformulating the convection term and the force term via the auxiliary variable, zero numerical flux and zero numerical source term are enforced at the hydrostatic equilibrium state instead of the balance between hydrostatic pressure and external force. Several problems are tested to demonstrate the accuracy and the stability of the new scheme. The results confirm that, the new scheme can preserve the exact hydrostatic solution. The small perturbation riding on hydrostatic equilibria can be calculated accurately. More importantly, the new scheme is capable of simulating the process of converging towards hydrostatic equilibria from a highly unbalanced initial condition. The ultimate state of zero velocity and constant temperature is achieved up to machine accuracy. As demonstrated by the numerical experiments, the current scheme is very suitable for small amplitude perturbation and long time running under gravitational potential.     

One-Dimensional Blood Flow with Discontinuous Properties and Transport: Mathematical Analysis and Numerical Schemes

In this paper we consider the one-dimensional blood flow model with discontinuous mechanical and geometrical properties, as well as passive scalar transport, proposed in [E.F. Toro and A. Siviglia. Flow in collapsible tubes with discontinuous mechanical properties: mathematical model and exact solutions. Communications in Computational Physics. 13(2), 361-385, 2013], completing the mathematical analysis by providing new propositions and new proofs of relations valid across different waves. Next we consider a first order DOT Riemann solver, proposing an integration path that incorporates the passive scalar and proving the well-balanced properties of the resulting numerical scheme for stationary solutions. Finally we describe a novel and simple well-balanced, second order, non-linear numerical scheme to solve the equations under study; by using suitable test problems for which exact solutions are available, we assess the well-balanced properties of the scheme, its capacity to provide accurate solutions in challenging flow conditions and its accuracy.     

A Computational Model for Simulation of Shallow Water Waves by Elastic Deformations in the Topography

We propose a coupled model to simulate shallow water waves induced by elastic deformations in the bed topography. The governing equations consist of the depth-averaged shallow water equations including friction terms for the water free-surface and the well-known second-order elastostatics formulation for the bed deformation. The perturbation on the free-surface is assumed to be caused by a sudden change in the bottom beds. At the interface between the water flow and the bed topography, transfer conditions are implemented. Here, the hydrostatic pressure and friction forces are considered for the elastostatic equations whereas bathymetric forces are accounted for in the shallow water equations. The focus in the present study is on the development of a simple and accurate representation of the interaction between water waves and bed deformations in order to simulate practical shallow water flows without relying on complex partial differential equations with free boundary conditions. The effects of location and magnitude of the deformation on the flow fields and free-surface waves are investigated in details. Numerical simulations are carried out for several test examples on shallow water waves induced by sudden changes in the bed. The proposed computational model has been found to be feasible and satisfactory.     

Finite Volume Evolution Galerkin Methods for the Shallow Water Equations with Dry Beds

We present a new Finite Volume Evolution Galerkin (FVEG) scheme for the solution of the shallow water equations (SWE) with the bottom topography as a source term. Our new scheme will be based on the FVEG methods presented in (Noelle and Kraft, J. Comp. Phys., 221 (2007)), but adds the possibility to handle dry boundaries. The most important aspect is to preserve the positivity of the water height. We present a general approach to ensure this for arbitrary finite volume schemes. The main idea is to limit the outgoing fluxes of a cell whenever they would create negative water height. Physically, this corresponds to the absence of fluxes in the presence of vacuum. Well-balancing is then re-established by splitting gravitational and gravity driven parts of the flux. Moreover, a new entropy fix is introduced that improves the reproduction of sonic rarefaction waves.     

A Space-Time Conservative Method for Hyperbolic Systems with Stiff and Non-Stiff Source Terms

In this article we propose a higher-order space-time conservative method for hyperbolic systems with stiff and non-stiff source terms as well as relaxation systems. We call the scheme a slope propagation (SP) method. It is an extension of our scheme derived for homogeneous hyperbolic systems [1]. In the present inhomogeneous systems the relaxation time may vary from order of one to a very small value. These small values make the relaxation term stronger and highly stiff. In such situations underresolved numerical schemes may produce spurious numerical results. However, our present scheme has the capability to correctly capture the behavior of the physical phenomena with high order accuracy even if the initial layer and the small relaxation time are not numerically resolved. The scheme treats the space and time in a unified manner. The flow variables and their slopes are the basic unknowns in the scheme. The source term is treated by its volumetric integration over the space-time control volume and is a direct part of the overall space-time flux balance. We use two approaches for the slope calculations of the flow variables, the first one results directly from the flux balance over the control volumes, while in the second one we use a finite difference approach. The main features of the scheme are its simplicity, its Jacobian-free and Riemann solver-free recipe, as well as its efficiency and high order accuracy. In particular we show that the scheme has a discrete analog of the continuous asymptotic limit. We have implemented our scheme for various test models available in the literature such as the Broadwell model, the extended thermodynamics equations, the shallow water equations, traffic flow and the Euler equations with heat transfer. The numerical results validate the accuracy, versatility and robustness of the present scheme.     

GPU Accelerated Discontinuous Galerkin Methods for Shallow Water Equations

We discuss the development, verification, and performance of a GPU accelerated discontinuous Galerkin method for the solutions of two dimensional nonlinear shallow water equations. The shallow water equations are hyperbolic partial differential equations and are widely used in the simulation of tsunami wave propagations. Our algorithms are tailored to take advantage of the single instruction multiple data (SIMD) architecture of graphic processing units. The time integration is accelerated by local time stepping based on a multi-rate Adams-Bashforth scheme. A total variational bounded limiter is adopted for nonlinear stability of the numerical scheme. This limiter is coupled with a mass and momentum conserving positivity preserving limiter for the special treatment of a dry or partially wet element in the triangulation. Accuracy, robustness and performance are demonstrated with the aid of test cases. Furthermore, we developed a unified multi-threading model OCCA. The kernels expressed in OCCA model can be cross-compiled with multi-threading models OpenCL, CUDA, and OpenMP. We compare the performance of the OCCA kernels when cross-compiled with these models.     

Localized Exponential Time Differencing Method for Shallow Water Equations: Algorithms and Numerical Study

In this paper, we investigate the performance of the exponential time differencing (ETD) method applied to the rotating shallow water equations. Comparing with explicit time stepping of the same order accuracy in time, the ETD algorithms could reduce the computational time in many cases by allowing the use of large time step sizes while still maintaining numerical stability. To accelerate the ETD simulations, we propose a localized approach that synthesizes the ETD method and overlapping domain decomposition. By dividing the original problem into many subdomain problems of smaller sizes and solving them locally, the proposed approach could speed up the calculation of matrix exponential vector products. Several standard test cases for shallow water equations of one or multiple layers are considered. The results show great potential of the localized ETD method for high-performance computing because each subdomain problem can be naturally solved in parallel at every time step.     

Solving Two-Mode Shallow Water Equations Using Finite Volume Methods

In this paper, we develop and study numerical methods for the two-mode shallow water equations recently proposed in [S. STECHMANN, A. MAJDA, and B. KHOUIDER, Theor. Comput. Fluid Dynamics, 22 (2008), pp. 407–432]. Designing a reliable numerical method for this system is a challenging task due to its conditional hyperbolicity and the presence of nonconservative terms. We present several numerical approaches–two operator splitting methods (based on either Roe-type upwind or central-upwind scheme), a central-upwind scheme and a path-conservative central-upwind scheme– and test their performance in a number of numerical experiments. The obtained results demonstrate that a careful numerical treatment of nonconservative terms is crucial for designing a robust and highly accurate numerical method.     

A New Approach of High Order Well-Balanced Finite Volume WENO Schemes and Discontinuous Galerkin Methods for a Class of Hyperbolic Systems with Source Terms

Hyperbolic balance laws have steady state solutions in which the flux gradients are nonzero but are exactly balanced by the source terms. In our earlier work [31–33], we designed high order well-balanced schemes to a class of hyperbolic systems with separable source terms. In this paper, we present a different approach to the same purpose: designing high order well-balanced finite volume weighted essentially non-oscillatory (WENO) schemes and RungeKutta discontinuous Galerkin (RKDG) finite element methods. We make the observation that the traditional RKDG methods are capable of maintaining certain steady states exactly, if a small modification on either the initial condition or the flux is provided. The computational cost to obtain such a well balanced RKDG method is basically the same as the traditional RKDG method. The same idea can be applied to the finite volume WENO schemes. We will first describe the algorithms and prove the well balanced property for the shallow water equations, and then show that the result can be generalized to a class of other balance laws. We perform extensive one and two dimensional simulations to verify the properties of these schemes such as the exact preservation of the balance laws for certain steady state solutions, the non-oscillatory property for general solutions with discontinuities, and the genuine high order accuracy in smooth regions.     

A Novel Strong-Coupling Pseudo-Spectral Method for Numerical Studies of Two-Layer Turbulent Channel Flows

A new method is proposed to simulate a coupled air-water two-layer turbulent channel flow. A numerically effective dynamic viscosity is implemented to calculate the viscous momentum flux at the interface, leading to a strong-coupling scheme for the evolution of air and water motions. The direct numerical simulation results are compared with those in the literature obtained from a weak-coupling scheme. It is discovered that while the turbulence statistics of the air phase based on the strong- and weak-coupling schemes are close to each other, the results on the water side are influenced by the coupling approach, especially near the water surface.     

High Order Finite Difference Hermite WENO Fixed-Point Fast Sweeping Method for Static Hamilton-Jacobi Equations

In this paper, we combine the nonlinear HWENO reconstruction in [43] and the fixed-point iteration with Gauss-Seidel fast sweeping strategy, to solve the static Hamilton-Jacobi equations in a novel HWENO framework recently developed in [22]. The proposed HWENO frameworks enjoys several advantages. First, compared with the traditional HWENO framework, the proposed methods do not need to introduce additional auxiliary equations to update the derivatives of the unknown function $\phi$. They are now computed from the current value of $\phi$ and the previous spatial derivatives of $\phi$. This approach saves the computational storage and CPU time, which greatly improves the computational efficiency of the traditional HWENO scheme. In addition, compared with the traditional WENO method, reconstruction stencil of the HWENO methods becomes more compact, their boundary treatment is simpler, and the numerical errors are smaller on the same mesh. Second, the fixed-point fast sweeping method is used to update the numerical approximation. It is an explicit method and does not involve the inverse operation of nonlinear Hamiltonian, therefore any Hamilton-Jacobi equations with complex Hamiltonian can be solved easily. It also resolves some known issues, including that the iterative number is very sensitive to the parameter $ε$ used in the nonlinear weights, as observed in previous studies. Finally, to further reduce the computational cost, a hybrid strategy is also presented. Extensive numerical experiments are performed on two-dimensional problems, which demonstrate the good performance of the proposed fixed-point fast sweeping HWENO methods.     

Simulation of Three-Dimensional Free-Surface Flows Using Two-Dimensional Multilayer Shallow Water Equations

We present an efficient and conservative Eulerian-Lagrangian method for solving two-dimensional hydrostatic multilayer shallow water flows with mass exchange between the vertical layers. The method consists of a projection finite volume method for the Eulerian stage and a method of characteristics to approximate the numerical fluxes for the Lagrangian stage. The proposed method is simple to implement, satisfies the conservation property and it can be used for multilayer shallow water equations on non-flat bathymetry including eddy viscosity and Coriolis forces. It offers a novel method of calculating stratified vertical velocities without the use of the Navier-Stokes equations. Numerical results are presented for several examples and the obtained results for a free-surface flow problem are in close agreement with the analytical solutions. We also test the performance of the proposed method for a test example of wind-driven flows with recirculation     

On Universal Osher-Type Schemes for General Nonlinear Hyperbolic Conservation Laws

This paper is concerned with a new version of the Osher-Solomon Riemann solver and is based on a numerical integration of the path-dependent dissipation matrix. The resulting scheme is much simpler than the original one and is applicable to general hyperbolic conservation laws, while retaining the attractive features of the original solver: the method is entropy-satisfying, differentiable and complete in the sense that it attributes a different numerical viscosity to each characteristic field, in particular to the intermediate ones, since the full eigenstructure of the underlying hyperbolic system is used. To illustrate the potential of the proposed scheme we show applications to the following hyperbolic conservation laws: Euler equations of compressible gasdynamics with ideal gas and real gas equation of state, classical and relativistic MHD equations as well as the equations of nonlinear elasticity. To the knowledge of the authors, apart from the Euler equations with ideal gas, an Osher-type scheme has never been devised before for any of these complicated PDE systems. Since our new general Riemann solver can be directly used as a building block of high order finite volume and discontinuous Galerkin schemes we also show the extension to higher order of accuracy and multiple space dimensions in the new framework of P_NP_M schemes on unstructured meshes recently proposed in [9].     

A One-Dimensional Second-Order Cell-Centered Lagrangian Scheme Satisfying the Entropy Condition

The numerical solutions of gas dynamics equations have to be consistent with the second law of thermodynamics, which is termed entropy condition. However, most cell-centered Lagrangian (CL) schemes do not satisfy the entropy condition. Until 2020, for one-dimensional gas dynamics equations, the first-order CL scheme with the hybridized flux developed by combining the acoustic approximate (AA) flux and the entropy conservative (EC) flux developed by Maire et al. was used. This hybridized CL scheme satisfies the entropy condition; however, it is under-entropic in the part zones of rarefaction waves. Moreover, the EC flux may result in nonphysical numerical oscillations in simulating strong rarefaction waves. Another disadvantage of this scheme is that it is of only first-order accuracy. In this paper, we firstly construct a modified entropy conservative (MEC) flux which can damp effectively numerical oscillations in simulating strong rarefaction waves. Then we design a new hybridized CL scheme satisfying the entropy condition for one-dimensional complex flows. This new hybridized CL scheme is a combination of the AA flux and the MEC flux.In order to prevent the specific entropy of the hybridized CL scheme from being under-entropic, we propose using the third-order TVD-type Runge-Kutta time discretization method. Based on the new hybridized flux, we develop the second-order CL scheme that satisfies the entropy condition.Finally, the characteristics of our new CL scheme using the improved hybridized flux are demonstrated through several numerical examples.     

Boundary Conditions for Limited Area Models Based on the Shallow Water Equations

A new set of boundary conditions has been derived by rigorous methods for the shallow water equations in a limited domain. The aim of this article is to present these boundary conditions and to report on numerical simulations which have been performed using these boundary conditions. The new boundary conditions which are mildly dissipative let the waves move freely inside and outside the domain. The problems considered include a one-dimensional shallow water system with two layers of fluids and a two-dimensional inviscid shallow water system in a rectangle.     

Consistent Forcing Scheme in the Simplified Lattice Boltzmann Method for Incompressible Flows

Considering the fact that the lattice discrete effects are neglected while introducing a body force into the simplified lattice Boltzmann method (SLBM), we propose a consistent forcing scheme in SLBM for incompressible flows with external forces. The lattice discrete effects are considered at the level of distribution functions in the present forcing scheme. Consequently, it is more accurate compared with the original forcing scheme used in SLBM. Through Taylor series expansion and Chapman-Enskog (CE) expansion analysis, the present forcing scheme can be proven to recover the macroscopic Navier-Stokes (N-S) equations. Then, the macroscopic equations are resolved through a fractional step technique. Furthermore, the material derivative term is discretized by the central difference method. To verify the results of the present scheme, we simulate with multiple forms of external force interactions including the space- and time-dependent body forces. Hence, the present forcing scheme overcomes the disadvantages of the original forcing scheme and the body force can be accurately imposed in the present scheme even when a coarse mesh is applied while the original scheme fails. Excellent agreements between the analytical solutions and our numerical results can be observed.     

Discrete-Velocity Vector-BGK Models Based Numerical Methods for the Incompressible Navier-Stokes Equations

In this paper, we propose a class of numerical methods based on discrete-velocity vector-BGK models for the incompressible Navier-Stokes equations. By analyzing a splitting method with Maxwell iteration, we show that the usual lattice Boltzmann discretization of the vector-BGK models provides a good numerical scheme. Moreover, we establish the stability of the numerical scheme. The stability and second-order accuracy of the scheme are validated through numerical simulations of the two-dimensional Taylor-Green vortex flows. Further numerical tests are conducted to exhibit some potential advantages of the vector-BGK models, which can be regarded as competitive alternatives of the scalar-BGK models.     

A Conservative Lagrangian Scheme for Solving Compressible Fluid Flows with Multiple Internal Energy Equations

Lagrangian methods are widely used in many fields for multi-material compressible flow simulations such as in astrophysics and inertial confinement fusion (ICF), due to their distinguished advantage in capturing material interfaces automatically. In some of these applications, multiple internal energy equations such as those for electron, ion and radiation are involved. In the past decades, several staggered-grid based Lagrangian schemes have been developed which are designed to solve the internal energy equation directly. These schemes can be easily extended to solve problems with multiple internal energy equations. However, such schemes are typically not conservative for the total energy. Recently, significant progress has been made in developing cell-centered Lagrangian schemes which have several good properties such as conservation for all the conserved variables and easiness for remapping. However, these schemes are commonly designed to solve the Euler equations in the form of the total energy, therefore they cannot be directly applied to the solution of either the single internal energy equation or the multiple internal energy equations without significant modifications. Such modifications, if not designed carefully, may lead to the loss of some of the nice properties of the original schemes such as conservation of the total energy. In this paper, we establish an equivalency relationship between the cell-centered discretizations of the Euler equations in the forms of the total energy and of the internal energy. By a carefully designed modification in the implementation, the cell-centered Lagrangian scheme can be used to solve the compressible fluid flow with one or multiple internal energy equations and meanwhile it does not lose its total energy conservation property. An advantage of this approach is that it can be easily applied to many existing large application codes which are based on the framework of solving multiple internal energy equations. Several two dimensional numerical examples for both Euler equations and three-temperature hydrodynamic equations in cylindrical coordinates are presented to demonstrate the performance of the scheme in terms of symmetry preserving, accuracy and non-oscillatory performance.