Simple and High-Accurate Schemes for Hyperbolic Conservation Laws

The paper constructs a class of simple high-accurate schemes (SHA schemes) with third order approximation accuracy in both space and time to solve linear hyperbolic equations, using linear data reconstruction and Lax-Wendroff scheme. The schemes can bemade even fourth order accurate with special choice of parameter. In order to avoid spurious oscillations in the vicinity of strong gradients, we make the SHA schemes total variation diminishing ones (TVD schemes for short) by setting flux limiter in their numerical fluxes and then extend these schemes to solve nonlinear Burgers’ equation and Euler equations.The numerical examples show that these schemes give high order of accuracy and high resolution results.The advantages of these schemes are their simplicity and high order of accuracy.


Introduction
In designing numerical schemes of very high order of accuracy for solving hyperbolic conservation laws one faces at least three major difficulties.One of them concerns the preservation of high accuracy in both space and time for multidimensional problems.Another one concerns conservation; this is mandatory in the presence of shock waves.The other very important issue relates to the generation of spurious oscillations in the vicinity of strong gradients.According to Godunov's theorem [1] these are unavoidable by linear schemes of accuracy greater than one.Each one of these difficulties is in itself not easy to resolve; the simultaneous resolution of all three difficulties above represents a formidable task in the numerical analysis of hyperbolic conservation laws.At present, there are various approaches for constructing numerical schemes that attempt to overcome the above difficulties.State-of-the-art very high order methods (at least third order) for hyperbolic conservation laws include the class of ENO/WENO schemes [2][3][4], Spectral Method [5], the class of Compact Difference Methods [6], Discontinuous Galerkin Finite Element Methods [7], and ADER Methods [8].These schemes meet the requirement of very high order spatial accuracy; matching time accuracy to space accuracy, however, remains an issue in all of the above approaches.
As regards the ENO/WENO/MPWENO approach, the most accurate scheme reported so far uses 9th order spatial discretisation with Runge-Kutta methods for time evolution.To avoid spurious oscillations these Runge-Kutta methods must be TVD.This leads to accuracy barriers: the accuracy of such methods cannot be higher than fifth.Moreover, fourth and fifth order methods are quite complicated and have reduced stability range.In most practical implementation reported, when the solution is not smooth, a third order TVD Runge-Kutta method has been used, for example, [2].Although increased order of spatial discretization improves accuracy for some problems, such schemes converge with third order only when the mesh is refined.Therefore when numerical simulation of compressible turbulence and wave propagation problems involves long-time evolution, it would be beneficial to use schemes which are simple and efficient and converge with higher order both in time and space.
How to construct simple and high-accurate schemes in both space and time by simple approach is an interesting work.There are various ways to construct high order numerical schemes for conservation laws.One of them is the MUSCL (Monotone Upstream-Centred Scheme for Conservation Laws) [9][10][11] approach which is a very popular and simple one to construct explicit, fully discrete schemes that lies on the solution of Riemann problem and could be extended to nonlinear systems.The MUSCL-type methods based on piecewise linear, local reconstruction of data are 2nd order of accuracy and their available highest order is 3rd for linear advection equation.In this paper, we generalize the MUSCL approach by modifying the piecewise linear reconstruction with a gradual slope, solving a Riemann problem at each interface and using Lax-Wendroff approach to compute the numerical flux.The new schemes constructed can reach 3rd order of accuracy in time and space, and its available highest order of accuracy is 4th.The flux limiter approach is used to obtain TVD version of the new schemes to avoid unphysical oscillations in the vicinity of large gradients.
The rest of the paper is organized as follows.In Section 2 we describe the construction of the class of schemes SHA for linear advection equation and prove its convergent order; next we set the limiter to the fluxes of these schemes to make them TVD ones in Section 3 and carry out some numerical examples in Section 4 and then extend them to solve nonlinear Burgers' equation in Section 5 and carry out some numerical examples in Section 6 and finally apply the schemes to solve one-dimensional Euler equations in Section 7. The conclusions are drawn in Section 8.

The Construction of Simple High-Accurate Schemes (SHA Schemes) for Advection Equations
In this section we consider the Initial Boundary Value Problem (IBVP) for the linear advection equation in the domain [0, ]×[0, ] on the − plane.This consists of a PDE together with initial condition (IC) and boundary conditions (BCs); namely, with suitable boundary conditions.A mesh for discretizing the domain is a regular grid of dimension Δ (spacing of grid points in the -direction) by Δ (spacing in the -direction).We suppose that [0, ] is discretized by  + 1 equally spaced grid points.Then The mesh points of the mesh are positioned at (Δ, Δ) on the  −  plane, with  = 0, 1, . . .,  and  = 0, 1, . ... We will use the notation   = Δ,   = Δ.The discrete values of the function (, ) at (Δ, Δ) will be denoted by We will also use the symbol    to denote an approximation to the exact mesh value (Δ, Δ).The distinction will be made at appropriate places.
We could construct a class of simple high-accurate schemes (SHA schemes) for advection equation as follows.
Step 1. Data reconstruction with boundary extrapolated values We choose    =    − (1/3)Δ  ,    =    + (1/3)Δ  , where Step 2. Consider evolution of    ,    by a time (1/2)Δ according to Step 3. Solution of piecewise constant data Riemann problem This step solves the conventional Riemann problem at the interface  +1/2 with data the solution of which for linear advection equation ( 1) is where  = (Δ/Δ).Then applying the solution (9) to the Lax-Wendroff flux for linear advection equation we get a new flux After some algebra operations,  SHA +1/2 becomes Journal of Applied Mathematics 3 Also, we can substitute (9), the solution of (7), into the following type of  −  flux; or two-step Richtmyer type flux of  −  method to get the flux as (12), where    = (   ),   +1 = (  +1 ).
Conservative Scheme.Substitution of  SHA +1/2 as numerical flux into conservative method gives a new scheme: where As for linear advection equation () = , we have The following introduces a theorem on how to verify the accuracy of a linear scheme: Theorem 1 (Roe [12]).A scheme of the form (17) is ℎ ( ≥ 0) order accurate in space and time if and only if According to Roe's theorem, it is easy to verify that scheme (16) is third order accurate in space and time for any value  ∈ [−1, 1] and the scheme is fourth order accurate when  = /2.

TVD Version of SHA Scheme
According to Godunov's Theorem [1] there are no monotone, linear schemes (17) of second or higher order of accuracy.This means that any linear scheme of second or higher is nonmonotone and there will be spurious oscillation in the vicinity of high gradients.The scheme must be modified in order to avoid spurious oscillation near the high gradients.
Harten [13] proposed the TVD property of conservation schemes.The numerical solutions will be free from spurious oscillations if the scheme has the TVD property.Now we will modify the SHA method above to a scheme with the TVD property according to Harten's Theorem.
Considering the class of nonlinear schemes where Δ +1/2 =  +1 −  and the coefficients  −1/2 ,  +1/2 are in general assumed to be functions of the data {   }, Harten had proved the following important result.
Theorem 2 (Harten [13]).For any scheme of the form (19) to solve (1), a sufficient condition for the scheme to be TVD is that the coefficients satisfy The conservative scheme SHA can be written as follows when applied to linear advection equation: where Thus we obtain Let or Then we could modify (23) as or where  is a function Obviously ( −1/2 ) (or ( +1/2 )) ∈ [0, 1]; schemes (27a) and (27b) modified above are TVD according to Harten's result.In addition formula (27a) can be rewritten into Formula (27b) also can be done similarly.So it can be seen that the conservativity of the TVD schemes (27a) and (27b) is not destroyed.

Numerical Experiments for
Linear Advection Equation 4.1.Accuracy Test for SHA Scheme.In the subsection we first test the order of accuracy of scheme SHA (16).Let us consider the linear equation subject to periodic initial data: This problem admits the global smooth solution that was computed at time  = 1 with the varying number of grid points .In Table 1 where () = .
In this experiment the Courant number  is  = 0.9. Figure 1 shows the excellent behavior of scheme (16) after long time.Figures 3 and 5 show the results computed by scheme (27a) and scheme (27b), respectively, at  = 9 units and reveal that schemes ((27a)-(27b)) also have high approximate accuracy to smooth initial problem.

Numerical Experiment II.
We solve the same problem as in experiment I, except that we replace the smooth initial condition with a discontinuous initial condition: where () = .
In experiment II the  is still 0.9 and Δ = 1/80.Figure 2 shows the result computed by the 4th order scheme SHA (16) ( = /2) at  = 2 units and reveals the 4th order scheme SHA without TVD products spurious oscillation in the vicinity of discontinuity. Figure 4 shows the result computed by the scheme SHA with TVD (27a) at  = 2 units and Figure 6 shows the result computed by scheme (27b) at  = 2 units.From the results shown by Figures 4 and 6 we see that schemes (27a) and (27b) both control the spurious oscillation, but the comparison of Figures 4 and 6 reveals that scheme (27b) has better resolution to discontinuous initial problem.The period boundary conditions are used in both experiments I and II.

The Application of SHA Scheme to Burgers' Equation
We can see the SHA schemes perform well in solving linear advection equation from the numerical results above.It can reach fourth order accuracy in space and time where the solution is smooth enough and the implementation of the scheme is so easy that the method can be extended to nonlinear conservative laws conveniently.In the following we will extend the schemes to inviscid Burgers' equation 2 ) with suitable boundary conditions.which is /.The full discrete conservative formula for the above equation is In the application of SHA scheme to the inviscid Burgers' equation we use the following transmissive boundary conditions: to provide numerical fluxes  1/2 ,  +1/2 .Denoting by   max the maximum wave speed throughout the domain at time level  we define the maximum Courant number   : where Δ is such that The time step Δ follows as where The calculation of the max wave speed value in (39) needs to use all of the speed values in the cells, including the boundary cells.

SHA Schemes for Solving Nonlinear Scalar Equation.
The high order scheme is constructed as follows.
Step 1. Data reconstruction Δ  is the slope on cell  for   (), and Then we get the values on the interface of cells: Step (45) Step 3. Solution of piecewise constant data Riemann problem Step 4. Compute the integral Then the numerical flux is computed by  SHA +1/2 = ( +1/2 +1/2 ).
The  +1/2 +1/2 is required for computing the numerical flux  SHA +1/2 ; that means the integral of the Riemann solution needs to be calculated.The following gives the exact integral with exact Riemann solution of inviscid Burgers' equation.There are two cases.

TVD Version of SHA for Burgers' Equation.
In this section we will construct the TVD version of SHA (52) by the combination of a low order monotone flux  LO +1/2 and a high order flux  HI +1/2 as where  +1/2 is a flux limiter.Here a specific Flux Limiter Centred scheme is constructed by taking where  force +1/2 is a First Order Centred scheme and is given as follows for nonlinear systems.
A classical scheme of first order accuracy to solve hyperbolic conservation laws is that of Lax-Friedrichs, whose numerical flux at the interface of two states   ,   is A second order accurate scheme Richtmyer scheme is given by computing a numerical flux where the intermediate state Hence  force +1/2 is given as As we know in [14], the relationship between conventional upwind flux limiters () and centred flux limiters () is with where   is the maximum Courant number.By choosing the upwind flux limiters SUPERBEE [15], VANLEER [16], VANALBADA [17], and MINBEE [15] we can obtain corresponding centred flux limiters.Here we recommend the use of the SUPERBEE flux limiter  sb to obtain a centred flux limiter: Set A single flux limiter is computed as Here  sb is the limiter function (61), and we will apply  LR to the flux in (53).As we see, the method can be easily extended to a system of conserved variables; we just need to apply (62)-(63) to every component   and obtain a corresponding limiter  LR  ,  = 1, . . ., .Then, one can select the final limiter as In next section, we present some numerical results for above methods.

Numerical Experiment for Inviscid Burgers' Equation
Apply the SHA scheme (52) to the inviscid Burgers' equation as follows: 2 ) the exact solution of which is  = sin ( − ).Firstly, we test the accuracy of the scheme when the solution is still continuous. 1 and  ∞ errors are listed in  Table 2 using SHA scheme (52) at  = 0.1 units.We have the order of accuracy, 2, in both  1 and  ∞ norms.
Next we compute the shock by SHA scheme (52) with  = 0.8 and Δ = 0.01.The comparison of exact solution and numerical solution for initial value problem (65) is given in Figure 9. Then we use the SHA scheme (53) with the limiter  LR () in (63) to solve the initial value problem (65).The comparison is shown in Figure 10.From Figures 9 and 10 we can see that SHA scheme gives satisfactory solution when shock develops at  = 0.4.
Considering the discontinuous initial value problem for Burgers' equation 2 ) The comparison of exact solution and numerical solution for initial value problem (66) computed by SHA scheme (52) is given in Figure 11; then we use the SHA scheme (53) with TVD to solve the initial value problem (66); the comparison is shown in Figure 12.From Figures 11 and 12 we see that SHA scheme (52) cannot control the spurious oscillation, but SHA scheme (53) not only can control the spurious oscillation but also has better resolution to discontinuity.

The Application of SHA Scheme to Euler Equations
We have presented every procedure of SHA scheme applying to Burgers' equation; now we will turn our focus on the application of SHA scheme to the following Euler equations: where and  is conservative variable,  is vector of convective flux,  is specific internal energy,  is total energy, and  is ratio of specific heats, the value of which is 1.4 for ideal gas.
The process of applying SHA to Euler equations is similar to that of Burgers' equation.The time step is determined by where   max is the maximum wave speed at the time   .A practical choice of   max is where    is the sound speed, for ideal gas    = √   /   .The data reconstruction and evolution part is the same as Burgers' equation.When we obtain the    and    at each interface, there are two ways to construct the flux  SHA +1/2 for each interface.
The Riemann solver free way is to use two-step Richtmyer version of Lax-Wendroff method to obtain the SHA numerical flux: where ).Then (53) can be used to construct TVD scheme for the Euler equations which is similar to the procedure for Burgers' equation.
This is a very common test for numerical scheme and its solution consists of a left rarefaction, a contact, and a right  This is a tough test for numerical schemes.The Courant number in this test is still 0.9 and Δ = 0.005.The numerical results computed by the five schemes are shown in Figures 16, 17, and 18 comparing with the exact solution.We can see that the numerical results from both SHA schemes are good, the results of SHA 2 are comparable with results obtained by WENO or DG, but SHA 2 is more simple and requires less computational cost.These show the effectiveness of our method.

Conclusion
This paper introduces a class of SHA schemes with 3rd order approximation accuracy only for linear advection equation by linear data reconstruction, solving Riemann problem and using Lax-Wendroff flux, and its available highest order is 4th.The general flux limiter approach is used to make the SHA schemes TVD in order to avoid spurious oscillations in the vicinity of strong gradients for nonlinear equations and systems.The numerical results show that these schemes can solve hyperbolic conservation laws with high order of accuracy, high resolution, and less computational cost.The next work is to extend SHA schemes to high-dimensional hyperbolic equations and apply them to compute some practical flows, for example, around an airfoil and wing.

Figure 1 : 2 Figure 2 :
Figure 1: SHA method applied to linear advection equation with a continuous initial value; numerical (symbols) and exact (line) solutions are compared at  = 9 units.

Figure 3 : 2 Figure 4 :
Figure3: SHA method with TVD (27a) applied to linear advection equation with a continuous initial value; numerical (symbols) and exact (line) solutions are compared at  = 9 units.

Figure 9 :
Figure 9: SHA method applied to inviscid Burgers' equation for continuous initial value; numerical (symbols) and exact (line) solutions are compared at  = 0.4 units.

Figure 10 :
Figure 10: SHA method with TVD applied to inviscid Burgers' equation for continuous initial value; numerical (symbols) and exact (line) solutions are compared at  = 0.4 units.

Figure 11 :
Figure 11: SHA method applied to inviscid Burgers' equation for discontinuous initial value; numerical (symbols) and exact (line) solutions are compared at  = 0.15 units.

Figure 15 :Figure 16 :
Figure 15: Pressure comparison of schemes at time 0.2 units for Sod's tube problem.

Table 1 :
Initial value problem (30) and  1 and  ∞ norms of the errors.
the accuracy of the scheme is shown.The table lists the errors at  = 1 unit using the fourth order scheme SHA (16) ( = /2); notice that we have the full order of accuracy, 4, in both  1 and  ∞ norms.
(1).Numerical Experiment I.Next we will use the fourth order scheme SHA (16) ( = /2) and TVD schemes (27a) and (27b) to compute the continuous initial value test problem (31) and discontinuous initial value test problem (32) for the linear advection(1).Computed results are shown, respectively, in Figures1, 3, and 5 and Figures2, 4, and 6.The first test problem to be solved is

Table 2 :
Initial value problem (65) and  1 and  ∞ norms of the errors.
Figure12: SHA method with TVD applied to inviscid Burgers' equation for continuous initial value; numerical (symbols) and exact (line) solutions are compared at  = 0.15 units..The mesh used is Δ = 0.005 and the Courant number is 0.9.The comparisons between the numerical solution obtained by the five schemes and the exact solution are shown in Figures13, 14, and 15.We can see that both kinds of SHA schemes with TVD give satisfactory results.The numerical results from the SHA 1 without Riemann solver are similar to that of MUSCL-Hancock methods, and both methods give more dissipation than the other three methods, and SHA 2 Figure 14: Velocity comparison of schemes at time 0.2 units for Sod's tube problem.
x u shock