The anti‐dissipative,non‐monotone behavior of Petrov–Galerkin upwinding

The Petrov–Galerkin method has been developed with the primary goal of damping spurious oscillations near discontinuities in advection dominated flows. For time‐dependent problems, the typical Petrov–Galerkin method is based on the minimization of the dispersion error and the simultaneous selective addition of dissipation. This optimal design helps to dampen the oscillations prevalent near discontinuities in standard Bubnov–Galerkin solutions. However, it is demonstrated that when the Courant number is less than 1, the Petrov–Galerkin method actually amplifies undershoots at the base of discontinuities. This is shown in an heuristic manner, and is demonstrated with numerical experiments with the scalar advection and Richards' equations. A discussion of monotonicity preservation as a design criterion, as opposed to phase or amplitude error minimization, is also presented. The Petrov–Galerkin method is further linked to the high‐resolution, total variation diminishing (TVD) finite volume method in order to obtain a monotonicity preserving Petrov–Galerkin method.     

A multidimensional HLL‐Riemann solver for non‐linear hyperbolic systems

This article presents a numerical model that enables to solve on unstructured triangular meshes and with a high order of accuracy, Riemann problems that appear when solving hyperbolic systems. For this purpose, we use a MUSCL‐like procedure in a ‘cell‐vertex’ finite‐volume framework. In the first part of this procedure, we devise a four‐state bi‐dimensional HLL solver (HLL‐2D). This solver is based upon the Riemann problem generated at the barycenter of a triangular cell, from the surrounding cell‐averages. A new three‐wave model makes it possible to solve this problem, approximately. A first‐order version of the bi‐dimensional Riemann solver is then generated for discretizing the full compressible Euler equations. In the second part of the MUSCL procedure, we develop a polynomial reconstruction that uses all the surrounding numerical data of a given point, to give at best third‐order accuracy. The resulting over determined system is solved by using a least‐square methodology. To enforce monotonicity conditions into the polynomial interpolation, we use and adapt the monotonicity‐preserving limiter, initially devised by Barth (AIAA Paper 90‐0013, 1990). Numerical tests and comparisons with competing numerical methods enable to identify the salient features of the whole model. Copyright © 2010 John Wiley & Sons, Ltd.     

An efficient low‐dissipative WENO filter scheme

This paper presents a simple and efficient procedure developed for tracing discontinuities in flow fields. Numerical experiments are carried out to test the new sensor coupled with the associated nonlinear WENO dissipation filter developed to suppress the numerical dissipation. The tests show that, for a problem containing shocks and vortices, the implementation of the new sensor and an optimized WENO scheme can obtain a stable and highly resolved solution. The numerical experiments demonstrated that the new filter scheme performs efficiently both in parallel and serial running for one‐dimensional inviscid flow problems. Direct numerical simulation of a Mach 5 turbulent boundary layer over a flat plate was carried out to demonstrate the applicability of the scheme to the DNS practices. Copyright © 2011 John Wiley & Sons, Ltd.     

A high‐order alternating direction implicit method for the unsteady convection‐dominated diffusion problem

A high‐order alternating direction implicit (ADI) method for solving the unsteady convection‐dominated diffusion equation is developed. The fourth‐order Padé scheme is used for the discretization of the convection terms, while the second‐order Padé scheme is used for the diffusion terms. The Crank–Nicolson scheme and ADI factorization are applied for time integration. After ADI factorization, the two‐dimensional problem becomes a sequence of one‐dimensional problems. The solution procedure consists of multiple use of a one‐dimensional tridiagonal matrix algorithm that produces a computationally cost‐effective solver. Von Neumann stability analysis is performed to show that the method is unconditionally stable. An unsteady two‐dimensional problem concerning convection‐dominated propagation of a Gaussian pulse is studied to test its numerical accuracy and compare it to other high‐order ADI methods. The results show that the overall numerical accuracy can reach third or fourth order for the convection‐dominated diffusion equation depending on the magnitude of diffusivity, while the computational cost is much lower than other high‐order numerical methods. Copyright © 2011 John Wiley & Sons, Ltd.     

Development of a class of multiple time‐stepping schemes for convection–diffusion equations in two dimensions

In this paper we present a class of semi‐discretization finite difference schemes for solving the transient convection–diffusion equation in two dimensions. The distinct feature of these scheme developments is to transform the unsteady convection–diffusion (CD) equation to the inhomogeneous steady convection–diffusion‐reaction (CDR) equation after using different time‐stepping schemes for the time derivative term. For the sake of saving memory, the alternating direction implicit scheme of Peaceman and Rachford is employed so that all calculations can be carried out within the one‐dimensional framework. For the sake of increasing accuracy, the exact solution for the one‐dimensional CDR equation is employed in the development of each scheme. Therefore, the numerical error is attributed primarily to the temporal approximation for the one‐dimensional problem. Development of the proposed time‐stepping schemes is rooted in the Taylor series expansion. All higher‐order time derivatives are replaced with spatial derivatives through use of the model differential equation under investigation. Spatial derivatives with orders higher than two are not taken into account for retaining the linear production term in the convection–diffusion‐reaction differential system. The proposed schemes with second, third and fourth temporal accuracy orders have been theoretically explored by conducting Fourier and dispersion analyses and numerically validated by solving three test problems with analytic solutions. Copyright © 2006 John Wiley & Sons, Ltd.     

High resolution Godunov‐type schemes with small stencils

Higher‐order Godunov‐type schemes have to cope with the following two problems: (i) the increase in the size of the stencil that make the scheme computationally expensive, and (ii) the monotony‐preserving treatments (limiters) that must be implemented to avoid oscillations, leading to strong damping of the solution, in particular linear waves (e.g. acoustic waves). When too compressive, limiting procedures may also trigger the instability of oscillatory numerical solutions (e.g. in advection–dispersion phenomena) via the artificial amplification of the shorter modes. The present paper proposes a new approach to carry out the reconstruction. In this approach, the values of the flow variable at the edges of the computational cells are obtained directly from the reconstruction within these cells. This method is applied to the MUSCL and DPM schemes for the solution of the linear advection equation. The modified DPM scheme can capture contact discontinuities within one computational cell, even after millions of time steps at Courant numbers ranging from 1 to values as low as 10^‐4. Linear waves are subject to negligible damping. Application of the method to the DPM for one‐dimensional advection–dispersion problems shows that the numerical instability of oscillatory solutions caused by the over compressive, original DPM limiter is eliminated. One‐ and two‐dimensional shallow water simulations show an improvement over classical methods, in particular for two‐dimensional problems with strongly distorted meshes. The quality of the computational solution in the two‐dimensional case remains acceptable even for mesh aspect ratios Δx/Δy as large as 10. The method can be extend to the discretization of higher‐order PDEs, allowing third‐order space derivatives to be discretized using only two cells in space. Copyright © 2004 John Wiley & Sons, Ltd.     

The composite finite volume method on unstructured meshes for the two‐dimensional shallow water equations

A composite finite volume method (FVM) is developed on unstructured triangular meshes and tested for the two‐dimensional free‐surface flow equations. The methodology is based on the theory of the remainder effect of finite difference schemes and the property that the numerical dissipation and dispersion of the schemes are compensated by each other in a composite scheme. The composite FVM is formed by global composition of several Lax–Wendroff‐type steps followed by a diffusive Lax–Friedrich‐type step, which filters out the oscillations around shocks typical for the Lax–Wendroff scheme. To test the efficiency and reliability of the present method, five typical problems of discontinuous solutions of two‐dimensional shallow water are solved. The numerical results show that the proposed method, which needs no use of a limiter function, is easy to implement, is accurate, robust and is highly stable. Copyright © 2001 John Wiley & Sons, Ltd.     

Entropy dissipation scheme and minimums‐increase‐and‐maximums‐decrease slope limiter

In this paper, we continue to study the entropy dissipation scheme developed in former. We start with a numerical study of the scheme without the entropy dissipation term on the linear advection equation, which shows that the scheme is stable and numerical dissipation and numerical dispersion free for smooth solutions. However, the numerical results for discontinuous solutions show nonlinear instabilities near jump discontinuities. This is because the scheme enforces two related conservation properties in the computation. With this study, we design a so‐called ‘minimums‐increase‐and‐maximums‐decrease’ slope limiter in the reconstruction step of the scheme and delete the entropy dissipation in the linear fields and reduce the entropy dissipation terms in the nonlinear fields. Numerical experiments show improvements of the designed scheme compared with the results presented in former. However, the minimums‐increase‐and‐maximums‐decrease limiter is still not perfect yet, and better slope limiters are still sought. Copyright © 2011 John Wiley & Sons, Ltd.     

High‐order numerical simulations of flow‐induced noise

In this paper, the flow/acoustics splitting method for predicting flow‐generated noise is further developed by introducing high‐order finite difference schemes. The splitting method consists of dividing the acoustic problem into a viscous incompressible flow part and an inviscid acoustic part. The incompressible flow equations are solved by a second‐order finite volume code EllipSys2D/3D. The acoustic field is obtained by solving a set of acoustic perturbation equations forced by flow quantities. The incompressible pressure and velocity form the input to the acoustic equations. The present work is an extension of our acoustics solver, with the introduction of high‐order schemes for spatial discretization and a Runge–Kutta scheme for time integration. To achieve low dissipation and dispersion errors, either Dispersion‐Relation‐Preserving (DRP) schemes or optimized compact finite difference schemes are used for the spatial discretizations. Applications and validations of the new acoustics solver are presented for benchmark aeroacoustic problems and for flow over an NACA 0012 airfoil. Copyright © 2010 John Wiley & Sons, Ltd.     

Some applications of the concept of minimized integrated exponential error for low dispersion and low dissipation

Several techniques to optimize parameters that regulate dispersion and dissipation effects in finite difference schemes have been devised in our previous works. They all use the notion that dissipation neutralizes dispersion. These techniques are the minimized integrated square difference error (MISDE) and the minimized integrated exponential error for low dispersion and low dissipation (MIEELDLD). It is shown in this work based on several numerical schemes tested that the technique of MIEELDLD is more accurate than MISDE to optimize the parameters that regulate dispersion and dissipation effects with the aim of improving the shock‐capturing properties of numerical methods. First, we consider the family of third‐order schemes proposed by Takacs. We use the techniques MISDE and MIEELDLD to optimize two parameters, namely, the cfl number and another variable which also controls dispersion and dissipation. Second, these two techniques are used to optimize a numerical scheme proposed by Gadd. Moreover, we compute the optimal cfl for some multi‐level schemes in 1D. Numerical tests for some of these numerical schemes mentioned above are performed at different cfl numbers and it is shown that the results obtained are dependent on the cfl number chosen. The errors from the numerical results have been quantified into dispersion and dissipation using a technique devised by Takacs. Finally, we make use of a composite scheme made of corrected Lax–Friedrichs and the two‐step Lax–Friedrichs schemes like the CFLF4 scheme at its optimal cfl number, to solve some problems in 2D, namely: solid body rotation test, acoustics and the circular Riemann problem. Copyright © 2011 John Wiley & Sons, Ltd.     

Robustness of the hybrid DRP‐WENO scheme for shock flow computations

This paper describes a new variant of hybrid scheme that is constructed by a wave‐capturing scheme and a nonoscillatory scheme for flow computations in the presence of shocks. The improved fifth‐order upwind weighted essentially nonoscillatory scheme is chosen to be conjugated with the seven‐point dispersion‐relation‐preserving scheme by means of an adaptive switch function of grid‐point type. The new hybrid scheme can achieve a better resolution than the hybrid scheme which is based on the classical weighted essentially scheme. Ami Harten's multiresolution analysis algorithm is applied to density field for detecting discontinuities and setting point values of the switch function adaptively. Moreover, the tenth‐order central filter is applied in smooth part of the flow field for damping dispersion errors. This scheme can promote overall computational efficiency and yield oscillation‐free results in shock flows. The resolution properties and robustness of the new hybrid scheme are tested in both 1D and 2D linear and nonlinear cases. It performs well for computing flow problems with rich structures of weak/strong shocks and large/small vortices, such as the shock‐boundary layer interaction problem in a shock tube, which illustrates that it is very robust and accurate for direct numerical simulation of gas‐dynamics flows. Copyright © 2011 John Wiley & Sons, Ltd.     

High‐order compact numerical schemes for non‐hydrostatic free surface flows

The development of a numerical scheme for non‐hydrostatic free surface flows is described with the objective of improving the resolution characteristics of existing solution methods. The model uses a high‐order compact finite difference method for spatial discretization on a collocated grid and the standard, explicit, single step, four‐stage, fourth‐order Runge–Kutta method for temporal discretization. The Cartesian coordinate system was used. The model requires the solution of two Poisson equations at each time‐step and tridiagonal matrices for each derivative at each of the four stages in a time‐step. Third‐ and fourth‐order accurate boundaries for the flow variables have been developed including the top non‐hydrostatic pressure boundary. The results demonstrate that numerical dissipation which has been a problem with many similar models that are second‐order accurate is practically eliminated. A high accuracy is obtained for the flow variables including the non‐hydrostatic pressure. The accuracy of the model has been tested in numerical experiments. In all cases where analytical solutions are available, both phase errors and amplitude errors are very small. Copyright © 2006 John Wiley & Sons, Ltd.     

Augmented upwind numerical schemes for the groundwater transport advection‐dispersion equation with local operators

When solute transport is advection‐dominated, the advection‐dispersion equation approximates to a hyperbolic‐type partial differential equation, and finite difference and finite element numerical approximation methods become prone to artificial oscillations. The upwind scheme serves to correct these responses to produce a more realistic solution. The upwind scheme is reviewed and then applied to the advection‐dispersion equation with local operators for the first‐order upwinding numerical approximation scheme. The traditional explicit and implicit schemes, as well as the Crank‐Nicolson scheme, are developed and analyzed for numerical stability to form a comparison base. Two new numerical approximation schemes are then proposed, namely, upwind–Crank‐Nicolson scheme, where only for the advection term is applied, and weighted upwind‐downwind scheme. These newly developed schemes are analyzed for numerical stability and compared to the traditional schemes. It was found that an upwind–Crank‐Nicolson scheme is appropriate if the Crank‐Nicolson scheme is only applied to the advection term of the advection‐dispersion equation. Furthermore, the proposed explicit weighted upwind‐downwind finite difference numerical scheme is an improvement on the traditional explicit first‐order upwind scheme, whereas the implicit weighted first‐order upwind‐downwind finite difference numerical scheme is stable under all assumptions when the appropriate weighting factor (θ) is assigned.     

Dispersion and stability analyses of the linearized two‐dimensional shallow water equations in boundary‐fitted co‐ordinates

In the present investigation, a Fourier analysis is used to study the phase and group speeds of a linearized, two‐dimensional shallow water equations, in a non‐orthogonal boundary‐fitted co‐ordinate system. The phase and group speeds for the spatially discretized equations, using the second‐order scheme in an Arakawa C grid, are calculated for grids with varying degrees of non‐orthogonality and compared with those obtained from the continuous case. The spatially discrete system is seen to be slightly dispersive, with the degree of dispersivity increasing with an decrease in the grid non‐orthogonality angle or decrease in grid resolution and this is in agreement with the conclusions reached by Sankaranarayanan and Spaulding (J. Comput. Phys., 2003; 184 : 299–320). The stability condition for the non‐orthogonal case is satisfied even when the grid non‐orthogonality angle, is as low as 30° for the Crank Nicolson and three‐time level schemes. A two‐dimensional wave deformation analysis, based on complex propagation factor developed by Leendertse (Report RM‐5294‐PR, The Rand Corp., Santa Monica, CA, 1967), is used to estimate the amplitude and phase errors of the two‐time level Crank–Nicolson scheme. There is no dissipation in the amplitude of the solution. However, the phase error is found to increase, as the grid angle decreases for a constant Courant number, and increases as Courant number increases. Copyright © 2003 John Wiley & Sons, Ltd.     

The advantages of using high‐order finite differences schemes in laminar–turbulent transition studies

This paper presents various finite difference schemes and compare their ability to simulate instability waves in a given flow field. The governing equations for two‐dimensional, incompressible flows were solved in vorticity–velocity formulation. Four different space discretization schemes were tested, namely, a second‐order central differences, a fourth‐order central differences, a fourth‐order compact scheme and a sixth‐order compact scheme. A classic fourth‐order Runge–Kutta scheme was used in time. The influence of grid refinement in the streamwise and wall normal directions were evaluated. The results were compared with linear stability theory for the evolution of small‐amplitude Tollmien–Schlichting waves in a plane Poiseuille flow. Both the amplification rate and the wavenumber were considered as verification parameters, showing the degree of dissipation and dispersion introduced by the different numerical schemes. The results confirmed that high‐order schemes are necessary for studying hydrodynamic instability problems by direct numerical simulation. Copyright © 2005 John Wiley & Sons, Ltd.     

Two‐dimensional modeling for stability analysis of two‐phase stratified flow

The effect of wavelength and relative velocity on the disturbed interface of two‐phase stratified regime is modeled and discussed. To analyze the stability, a small perturbation is imposed on the interface. Growth or decline of the disturbed wave, relative velocity, and surface tension with respect to time will be discussed numerically. Newly developed scheme applied to a two‐dimensional flow field and the governing Navier–Stokes equations in laminar regime are solved. Finite volume method together with non‐staggered curvilinear grid is a very effective approach to capture interface shape with time. Because of the interface shape, for any time advancement, a new grid is performed separately on each stratified field, liquid, and gas regime. The results are compared with the analytical characteristics method and one‐dimensional modeling. This comparison shows that solving the momentum equation including viscosity term leads to physically more realistic results. In addition, the newly developed method is capable of predicting two‐phase stratified flow behavior more precisely than one‐dimensional modeling. It was perceived that the surface tension has an inevitable role in dissipation of interface instability and convergence of the two‐phase flow model. Copyright © 2009 John Wiley & Sons, Ltd.     

Finite‐difference 4th‐order compact scheme for the direct numerical simulation of instabilities of shear layers

An algorithm based on the 4th‐order finite‐difference compact scheme is developed and applied in the direct numerical simulations of instabilities of channel flow. The algorithm is illustrated in the context of stream function formulation that leads to field equation involving 4th‐order spatial derivatives. The finite‐difference discretization in the wall‐normal direction uses five arbitrarily spaced points. The discretization coefficients are determined numerically, providing a large degree of flexibility for grid selection. The Fourier expansions are used in the streamwise direction. A hybrid Runge–Kutta/Crank–Nicholson low‐storage scheme is applied for the time discretization. Accuracy tests demonstrate that the algorithm does deliver the 4th‐order accuracy. The algorithm has been used to simulate the natural instability processes in channel flow as well as processes occurring when the flow is spatially modulated using wall transpiration. Extensions to three‐dimensional situations are suggested. Copyright © 2005 John Wiley & Sons, Ltd.     

An implicit three‐dimensional fully non‐hydrostatic model for free‐surface flows

An implicit method is developed for solving the complete three‐dimensional (3D) Navier–Stokes equations. The algorithm is based upon a staggered finite difference Crank‐Nicholson scheme on a Cartesian grid. A new top‐layer pressure treatment and a partial cell bottom treatment are introduced so that the 3D model is fully non‐hydrostatic and is free of any hydrostatic assumption. A domain decomposition method is used to segregate the resulting 3D matrix system into a series of two‐dimensional vertical plane problems, for each of which a block tri‐diagonal system can be directly solved for the unknown horizontal velocity. Numerical tests including linear standing waves, nonlinear sloshing motions, and progressive wave interactions with uneven bottoms are performed. It is found that the model is capable to simulate accurately a range of free‐surface flow problems using a very small number of vertical layers (e.g. two–four layers). The developed model is second‐order accuracy in time and space and is unconditionally stable; and it can be effectively used to model 3D surface wave motions. Copyright © 2004 John Wiley & Sons, Ltd.     

Compact high‐resolution algorithms for time‐dependent advection on unstructured grids

A technique for constructing monotone, high resolution, multi‐dimensional upwind fluctuation distribution schemes for the scalar advection equation is presented. The method combines the second‐order Lax–Wendroff scheme with the upwind positive streamwise invariant (PSI) scheme via a fluctuation redistribution step, which ensures monotonicity (and which is a generalization of the flux‐corrected transport approach for fluctuation distribution schemes). Furthermore, the concept of a distribution point is introduced, which, when related to the equivalent equation for the scheme, leads to a ‘preferred direction’ for the limiting procedure, and hence to a new distribution of the fluctuation, which retains second‐order accuracy from the Lax–Wendroff scheme, even when the solution contains turning points. Experimental comparisons show that the new method compares favourably in terms of speed, accuracy and robustness with other, similar, techniques. Copyright © 2000 John Wiley & Sons, Ltd.     

A horizontally curvilinear non‐hydrostatic model for simulating nonlinear wave motion in curved boundaries

A horizontally curvilinear non‐hydrostatic free surface model that embeds the second‐order projection method, the so‐called θ scheme, in fractional time stepping is developed to simulate nonlinear wave motion in curved boundaries. The model solves the unsteady, Navier–Stokes equations in a three‐dimensional curvilinear domain by incorporating the kinematic free surface boundary condition with a top‐layer boundary condition, which has been developed to improve the numerical accuracy and efficiency of the non‐hydrostatic model in the standard staggered grid layout. The second‐order Adams–Bashforth scheme with the third‐order spatial upwind method is implemented in discretizing advection terms. Numerical accuracy in terms of nonlinear phase speed and amplitude is verified against the nonlinear Stokes wave theory over varying wave steepness in a two‐dimensional numerical wave tank. The model is then applied to investigate the nonlinear wave characteristics in the presence of dispersion caused by reflection and diffraction in a semicircular channel. The model results agree quantitatively with superimposed analytical solutions. Finally, the model is applied to simulate nonlinear wave run‐ups caused by wave‐body interaction around a bottom‐mounted cylinder. The numerical results exhibit good agreement with experimental data and the second‐order diffraction theory. Overall, it is shown that the developed model, with only three vertical layers, is capable of accurately simulating nonlinear waves interacting within curved boundaries. Copyright © 2011 John Wiley & Sons, Ltd.