Numerically induced oscillations in finite element approximations to the shallow water equations

Numerical noise has been a problem with finite element solutions to the shallow water equations. Two methods used to reduce the noise level are evaluated, and these results are compared with published results for equal-order interpolations. The two methods are mixed-interpolation (quadratic interpolation for velocity and linear interpolation for sea level) and a spectral form of the wave equation. Whereas mixed interpolation removes the troublesome sea level mode, it can still have considerable noise in velocity. The spectral wave equation is efficient and does not contain the spurious eigenmodes which contribute to high noise levels.     

Numerical simulation of liquefied fuel spills: I. Instantaneous release into a confined area

The shallow-water equations in radial symmetry are solved numerically to simulate the collapse of a cylindrical liquid column into an area surrounded by a concentric dike. The following three subcases of this problem are considered: a liquid column collapsing onto a layer of the same liquid, a liquid column collapsing onto a solid surface, and a column of lighter liquid collapsing onto a heavier liquid (i.e. liquefied natural gas (LNG) spilled onto water). The results for the three categories are compared and the differences and similarities between them are analysed.     

Solution of shallow water equations for complex flow domains via boundary-fitted co-ordinates

Many problems of applied oceanography and environmental science demand the solution of the momentum, mass and energy equations on physical domains having curving coastlines. Finite-difference calculations representing the boundary as a step function may give inaccurate results near the coastline where simulation results are of greatest interest for numerous applications. This suggests the use of methods which are capable of handling the problem of boundary curvature. This paper presents computational results for the shallow water equations on a circular ring of constant depth, employing the concept of boundary fitted grids (BFG) for an accurate representation of the boundary. All calculations are performed on a rectangle in the transformed plane using a mesh with square grid spacing. Comparisons of the simulations of transient normal mode oscillations and analytic solutions are shown, demonstrating that this technique yields accurate results in situations (provided that there is a reasonable choice of grid) involving a curved boundary. The software developed allows application to any two-dimensional area, regardless of the complexity of the geometry. Simulation runs were made with two co-ordinate systems. For the first system, the grid point distribution was obtained from polar co-ordinates. For the second one, grid point positions were calculated numerically, solving Poisson's equation. It was found that small variations in the metric coefficients do not deteriorate the accuracy of the simulation results. Moreover, comparisons of surface elevation and velocity components at grid points near the inner and outer radii obtained from an x?y Cartesian grid model with the BFG simulation were made. The former model produced inacccuracies at grid points near boundaries, and, owing to the large number of mesh points used to yield the necessary fine resolution, the computation time was found to be a factor of three higher.     

Numerical simulation of liquefied fuel spills: II. Instantaneous and continuous LNG spills on an unconfined water surface

A simple numerical model based on the shallow water equations in radial symmetry is used to simulate both instantaneous and continuous spills of liquefied natural gas (LNG) onto a water surface. Using the computed results, a study is made of the similarities and differences in the pool structure resulting from the two types of spills. For instantaneous spills a relation linear on a logarithmic plot is suggested between the maximum pool size and the spill volume. The effects of shear forces and surface cohesivity on the evolution of the spill are also examined.     

The effect of using non-orthogonal boundary-fitted grids for solving the shallow water equations

The purpose of this paper is to examine the effects of using non-orthogonal boundary-fitted grids for the numerical solution of the shallow water equations. Two geometries with well known analytical solutions are introduced in order to investigate the accuracy of the numerical solutions. The results verify that a reasonable departure from orthogonality can be allowed when the rate of change of cell areas is kept sufficiently small (i.e. it is not necessary to create a strictly orthogonal grid when the grid is sufficiently smooth).     

River and reservoir flow modelling using the transformed shallow water equations

This paper describes a versatile finite difference scheme for the solution of the two-dimensional shallow water equations on boundary-fitted non-orthogonal curvilinear meshes. It is believed that this is the first non-orthogonal shallow water equation model incorporating the advective acceleration terms to have been developed in the United Kingdom. The numerical scheme has been validated against the severe condition of jet-forced flow in a circular reservoir with vertical side walls, where reflections of the initial free surface waves pose major problems in achieving a stable solution. Furthermore, the validation exercises are designed to test the computer model for artificial diffusion, which may be a consequence of the numerical scheme adopted to stabilize the shallow water equations. The model is shown to be capable of simulating the flow conditions in an irregularly shaped domain typical of the geometries frequently encountered in civil engineering river basin management.     

Approximation of shallow water equations by Roe's Riemann solver

The inviscid shallow water equations provide a genuinely hyperbolic system and all the numerical tools that have been developed for a system of conservation laws can be applied to them. However, this system of equations presents some peculiarities that can be exploited when developing a numerical method based on Roe's Riemann solver and enhanced by a slope limiting of MUSCL type. In the present paper a TVD version of the Lax-Wendroff scheme is used and its performance is shown in 1D and 2D computations. Then two specific difficulties that arise in this context are investigated. The former is the capability of this class of schemes to handle geometric source terms that arise to model the bottom variation. The latter analysis pertains to situations in which strict hyperbolicity is lost, i.e. when two eigenvalues collapse into one.     

Bidirectional Whitham equations as models of waves on shallow water

Hammack & Segur (1978) conducted a series of surface water-wave experiments in which the evolution of long waves of depression was measured and studied. This present work compares time series from these experiments with predictions from numerical simulations of the KdV, Serre, and five unidirectional and bidirectional Whitham-type equations. These comparisons show that the most accurate predictions come from models that contain accurate reproductions of the Euler phase velocity, sufficient nonlinearity, and surface tension effects. The main goal of this paper is to determine how accurately the bidirectional Whitham equations can model data from real-world experiments of waves on shallow water. Most interestingly, the unidirectional Whitham equation including surface tension provides the most accurate predictions for these experiments. If the initial horizontal velocities are assumed to be zero (the velocities were not measured in the experiments), the three bidirectional Whitham systems examined herein provide approximations that are significantly more accurate than the KdV and Serre equations. However, they are not as accurate as predictions obtained from the unidirectional Whitham equation.     

On the effects of irregular boundaries in finite difference models

Propagation of periodic waves in the vicinity of irregular saw-tooth shaped boundary in finite difference models is investigated. The reflection of an incoming wave from a single saw-tooth boundary is found to be accompanied by a phase shift. It is shown that any wave mode propagating along such a boundary is trapped and decays in the direction normal to the boundary. A wave propagating along a channel with saw-tooth shaped lateral boundaries is influenced by the trapped waves, which leads to a reduction of the phase velocity. Phase velocities obtained from the present normal mode analysis are compared to velocities in numerical experiments. The agreement is excellent.     

A study of the motion in rotating containers using a boundary fitted co-ordinate system

Numerical solutions are often inaccurate because conventional co-ordinate systems do not represent the complex physical boundaries accurately. In the present work, the numerical solution of linear shallow water wave equations has been obtained by transforming the physical domain into a rectangular computational domain using elliptic differential operators. This work is part of a programme to develop three-dimensional body-fit grid systems for environmental flows. Solutions have been obtained for a cylindrical container and also a parabolic container. The initial conditions chosen are the ones for which analytical solutions exist. The numerical solutions compare well with analytical solutions.     

Godunov-type solution of the shallow water equations on adaptive unstructured triangular grids

A Godunov-type upwind finite volume solver of the non-linear shallow water equations is described. The shallow water equations are expressed in a hyperbolic conservation law formulation for application to cases where the bed topography is spatially variable. Inviscid fluxes at cell interfaces are computed using Roe's approximate Riemann solver. Second-order accurate spatial calculations of the fluxes are achieved by enhancing the polynomial approximation of the gradients of conserved variables within each cell. Numerical oscillations are curbed by means of a non-linear slope limiter. Time integration is second-order accurate and implicit. The numerical model is based on dynamically adaptive unstructured triangular grids. Test cases include an oblique hydraulic jump, jet-forced flow in a flat-bottomed circular reservoir, wind-induced circulation in a circular basin of non-uniform bed topography and the collapse of a circular dam. The model is found to give accurate results in comparison with published analytical and alternative numerical solutions. Dynamic grid adaptation and the use of a second-order implicit time integration scheme are found to enhance the computational efficiency of the model.     

Numerical solutions of the shallow water equations with discontinuous bed topography

A simple scheme is developed for treatment of vertical bed topography in shallow water flows. The effect of the vertical step on flows is modelled with the shallow water equations including local energy loss terms. The bed elevation is denoted with z_b^‐ for the left and z_b⁺ for the right values at each grid point, hence exactly representing a discontinuity in the bed topography. The surface gradient method (SGM) is generalized to reconstruct water depths at cell interfaces involving a vertical step so that the fluxes at the cell interfaces can accurately be calculated with a Riemann solver. The scheme is verified by predicting a surge crossing a step, a tidal flow over a step and dam‐break flows on wet/dry beds. The results have shown good agreements compared with analytical solutions and available experimental data. The scheme is efficient, robust, and may be used for practical flow calculations. Copyright © 2002 John Wiley & Sons, Ltd.     

Modified shallow water equations which admit the propagation of discontinuous waves over a dry bed

A method for modeling the propagation of discontinuous waves over a dry bed using the first approximation of shallow water theory is proposed. The method is based on a modified conservation law of total momentum that takes into account the concentrated momentum losses due to the formation of local turbulent vortex structures in the fluid surface layer at a discontinuous-wave front. A quantitative estimate of these losses is obtained by deriving the shallow water equations from the Navier-Stokes equations with allowance for viscosity, which has a rapidly increasing effect in the turbulent flow regions described by discontinuous waves. The stability of the discontinuous waves admitted by the modified system of conservation laws of shallow water theory is examined. As an example, a comparative analysis is performed of the solutions of the dam-break problem obtained for the classical and modified shallow water models. __________ Translated from Prikladnaya Mekhanika i Tekhnicheskaya Fizika, Vol. 48, No. 6, pp. 22–43, November–December, 2007     

Numerical solution of shallow water magnetohydrodynamic equations with non-flat bottom topography

The governing equations of shallow water magnetohydrodynamics describe the dynamics of a thin layer of nearly incompressible and electrically conducting fluids for which the evolution is nearly two-dimensional with magnetic equilibrium in the third direction. A high-resolution central-upwind scheme is applied to solve the model equations considering non-flat bottom topography. The suggested method is an upwind biased non-oscillatory finite volume scheme which doées not require a Riemann solver at each time step. To satisfy the divergence-free constraint, the projection method is used. Several case studies are carried out. For validation, a gas kinetic flux vector splitting scheme is also applied to the same model.     

Numerical viscosity reduction in the resolution of the shallow water equations with turbulent term

A first‐order finite volume model for the resolution of the 2D shallow water equations with turbulent term is presented. An upwind discretization of the equations that include the turbulent term is carried out. A method to reduce the excess of numerical viscosity (or diffusion) produced by the upwinding of the flux term is proposed. Two different discretizations of the turbulent term are compared, and results for uniform distributions of the viscosity are presented. Finally, two discretizations of the time derivative which are more efficient than Euler's are proposed and compared. Copyright © 2008 John Wiley & Sons, Ltd.     

A second-order radiation boundary condition for the shallow water wave equations on two-dimensional unstructured finite element grids

A second-order radiation boundary condition (RBC) is derived for 2D shallow water problems posed in ‘wave equation’ form and is implemented within the Galerkin finite element framework. The RBC is derived by matching the dispersion relation for the interior wave equation with an approximate solution to the exterior problem for outgoing waves. The matching is correct to second order, accounting for curvature of the wave front and the geometry. Implementation is achieved by using the RBC as an evolution equation for the normal gradient on the boundary, coupled through the natural boundary integral of the Galerkin interior problem. The formulation is easily implemented on non-straight, unstructured meshes of simple elements. Test cases show fidelity to solutions obtained on extended meshes and improvement relative to simpler first-order RBCs.     

The influence of normal flow boundary conditions on spurious modes in finite element solutions to the shallow water equations

Finite element solutions of the primitive equation (PE) form of the shallow water equations are notorious for the severe spurious 2Δx modes which appear. Wave equation (WE) solutions do not exhibit these numerical modes. In this paper we show that the severe spurious modes in PE solutions are strongly influenced by essential normal flow boundary conditions in the coupled continuity-momentum system of equations. This is demonstrated through numerical examples that avoid the use of essential normal flow boundary conditions either by specifying elevation values over the entire boundary or by implementing natural flow boundary conditions in the weak weighted residual form of the continuity equation. Results from a series of convergence tests show that PE solutions are of nearly the same quality as WE solutions when spurious modes are suppressed by alternative specification of the boundary conditions. Network intercomparisons indicate that varying nodal support does not excite spurious modes in a solution, although it does enhance the spurious modes introduced when an essential normal flow boundary condition is used. Dispersion analysis of discrete equations for interior and boundary nodes offers an explanation of the observed solution behaviour. For certain PE algorithms a mixed situation can arise where the boundary nodes exhibit a monotonic (noise-free) dispersion relationship and the interior nodes exhibit a folded (noisy) dispersion relationship. We have found that the mixed situation occurs when all boundary nodes are specified elevation nodes (which are enforced as essential conditions in the continuity equation) or when specified flow boundary nodes are treated as natural boundary conditions in the continuity equation. In either case the effect is to generate a solution that is essentially free of noise. Apparently, the monotonic dispersion behaviour at the boundaries suppresses the otherwise noisy behaviour caused by the folded dispersion relation on the interior.