An Improved Non-linear Weights for Seventh-Order WENO Scheme

An Improved Non-linear Weights for Seventh-Order WENO Scheme

Samala Rathan,  G Naga Raju
Department of Mathematics, Visvesvaraya National Institute of Technology, Nagpur, India
Email: rathan.maths@gmail.comEmail: gnagaraju@mth.vnit.ac.in
Abstract

In this article, the construction and implementation of a seventh order weighted essentially non-oscillatory scheme is reported for hyperbolic conservation laws. Local smoothness indicators are constructed based on -norm, where a higher order interpolation polynomial is used with each derivative being approximated to the fourth order of accuracy with respect to the evaluation point. The global smoothness indicator so constructed ensures the scheme achieves the desired order of accuracy. The scheme is reviewed in the presence of critical points and verified the numerical accuracy, convergence with the help of linear scalar test cases. Further, the scheme is implemented to non-linear scalar and system of equations in one and two dimensions. As the formulation is based on method of lines, to move forward in time linear strong-stability-preserving Runge-Kutta scheme for the linear problems and the fourth order nonlinear version of five stage strong stability preserving Runge-Kutta scheme for nonlinear problems is used.

Keywords: Hyperbolic conservation laws, non-linear weights, smoothness indicators, WENO scheme, Runge-Kutta schemes.
MSC: 65M20, 65N06, 41A10.

1 Introduction

The hyperbolic conservation laws arise in many applications such as in gas dynamics, magnetohydrodynamics (MHD) and shallow water flows. It is well known that even if the initial conditions are smooth, the hyperbolic conservation laws may develop discontinuities in its solution, such as shocks, contact discontinuities etc. Godunov [9] was first to propose a first order upwind scheme for the solution of these equations in the year which turned out to be a stepping stone for the development of various upwind schemes in the following years. In order to construct a higher order scheme Harten [12, 13] introduced the concept of Total Variation Diminishing (TVD), which says that the total variation to the approximation of numerical solution must be non-increasing with time. Later it was shown that the TVD schemes are having at most of first-order accuracy near smooth extrema [23] .

Harten et al. [14, 15, 16] derived the higher order schemes with the property of relaxing the TVD condition and allowing the occurrence of spurious oscillations in the numerical scheme but the Gibbs-like phenomena is essentially prevented, which is termed as Essentially Non-oscillatory (ENO) property and the schemes are known as ENO schemes. These are the first successful higher order schemes for the spatial discretization of the hyperbolic conservation laws, in a finite-volume formulation. The ENO schemes adopts a strategy of choosing the interpolation points over a stencil which avoids the induction of oscillations in the numerical solution through a smoothness indicator of a solution. And based on this idea the smoothest stencil is chosen from a set of candidate stencils. As a result, the ENO scheme obtains information from smooth regions and avoids spurious oscillations near discontinuities. Further, these schemes were studied in a finite-difference environment by Shu and Osher [31, 32].

The weighted ENO (WENO) schemes are set forth by Liu et al. [22], in a finite-volume frame of reference up to third-order of accuracy. Later, Jiang and Shu [19] have put forward these WENO schemes in a finite-difference setup to a higher order accuracy with the new smoothness indicators. These smoothness indicators are measured in the scaled norm, that is, they are the sum of the normalized squares, of all derivatives of the local interpolating polynomials. This scheme is referred as WENO-JS in the content to follow. For more details on ENO and WENO schemes, one can refer to the articles [29] and [30]. A very high order schemes are constructed in a similar manner of WENO-JS in [2], which we mention them here as WENO-BS schemes. Seventh order WENO-BS scheme is revised in [27, 28] and inspected the scheme in the presence of critical points.

Henrick et al. [17] examined that the actual convergence rate of the fifth-order WENO-JS scheme is less than the desired order, for the problems where the first and third order derivatives of the flux do not vanish simultaneously. In addition, it was ascertained that the convergence rate of the scheme is sensitive to the parameter employed in the evaluation of smoothness indicators to overcome from vanishing denominator. The authors revived the WENO-JS scheme by using a mapping function on the nonlinear weights such that the scheme, named as mapped WENO, satisfies the sufficient condition where WENO-JS fails and achieves an optimal order of convergence near simple smooth extrema. Subsequently, a very high order WENO schemes were developed based on the mapping function in [8].

Borges et al. [3] reviewed the fifth-order WENO schemes, entitled as WENO-Z scheme, by initiating a global smoothness indicator, which measures the smoothness of the larger stencil utilized in the construction of nonlinear weights. It was numerically validated that WENO-Z scheme is less dissipative than WENO-JS scheme and more efficient than mapped WENO scheme. WENO-Z scheme retained the convergence order as four at the first-order critical points, degrade to two when higher order critical points are encountered. These thoughts are extended by Castro et al. [4] to higher order schemes and produced a closed-form formula for the global smoothness indicators. The authors also assessed the dominance of the parameters and to retain the desired order of accuracy. The parameter is set up in the formulation of nonlinear weights to ascertain that these nonlinear weights converge to the ideal weights at a fast enough convergence rate.

The convergence analysis of WENO-JS scheme explored by Arandiga et al. [1] is based on the value of proposed that value is proportional to the square of mesh size , instead of a constant value so that the scheme achieves order of accuracy at smooth regions regardless of neighboring extrema, while this is of order when the function has a discontinuity in the stencil of points and is smooth in at least one of the point stencil. A question about the behavior of WENO-Z scheme when the value is taken in accordance with the value mentioned in [1], is examined by Don and Borges [5]. The authors made the accuracy analysis of the WENO-Z scheme and suggested a condition on the value of to achieve the full global-order of accuracy as similar to that of [1]. Further the authors have shown that the numerical oscillations can be attenuated by increasing the parameter value from to .

An alternate to the smoothness indicators of fifth order WENO-JS scheme were formulated by Fan et al. [6] with the help of Lagrange interpolation polynomials, accordingly a very high order schemes were derived by Fan in [7]. These schemes are based on the idea of constructing higher order global smoothness indicators, due to which less dissipation occurs in the solution near discontinuities. Very recently, another version of a smoothness indicators were proposed by Ha et al. [11], measured in -norm, hereafter referred as WENO-NS scheme. The authors introduced a higher-order approximation to the first derivative in the formation of local smoothness indicators which yields an improved behavior relative to other fifth-order WENO schemes. The global smoothness indicator for the WENO-NS scheme is preferred as an average of the two measurements, the smoothness information of the five point stencil and the middle three point stencil.

Kim et al. [20] perceived that, the three sub-stencils of the fifth-order WENO-NS scheme provides an unbalanced contribution to the flux at an evaluation point along the interface and an additional contribution term which measures the smoothness of the middle stencil in the formation of global smoothness indicator. The authors made a balanced tradeoff among the sub-stencils through a parameter and a global smoothness indicator is figured out which doesn’t depend on the smoothness information of the middle three point stencil anymore. These modifications lead to better results than the WENO-NS as well as other fifth-order WENO schemes, this scheme is pointed out as WENO-P scheme in the later part of this article.

A simple analysis verifies that the WENO-NS and WENO-P schemes attain the fifth-order accuracy at first order critical points but fails to achieve the accuracy at the points where second derivative vanishes. We have suggested a modified WENO-P scheme in [24] based on the idea of the linear combination of second-order derivatives, leading to a higher-order derivative information, is used in the construction of a global smoothness indicator. The modified smoothness indicator satisfies the sufficient condition, assert the requirement to achieve desired order of accuracy, even in the case of second order derivative vanishes.

In this article, a seventh order WENO scheme is derived in the lines of [11] and [24]. The smoothness indicators are obtained from the generalized undivided difference operator. Each of this operator is up to fourth-order of accuracy at the evaluation point, so the resulting scheme is seventh order accurate. Introduced parameters , to balance the tradeoff between the accuracies around the smooth to the discontinuous regions. The global smoothness indicator so earned satisfies the sufficient condition to get the optimal order of convergence rate, unvarying in the presence of critical points. Utilized strong stability preserving Runge-Kutta schemes introduced in [10] to advance the time. These are detailed out in the following sections, briefly they are:

Section , deals with the preliminaries about WENO reconstruction to the one-dimensional scalar conservation laws and section introduces the proposed WENO scheme where a new global measurement and local smoothness indicators are derived, which estimates the smoothness of a local solution in the construction of a seventh-order WENO scheme. Numerics for one-dimensional scalar test problems such as linear advection, Burger’s equation, the examples pertaining to the system of Euler equations such as shock tube problems, 1D shock-entropy wave interaction problem and 2D Riemann problem of gas dynamics are reported in Section , to demonstrate the advantages of the proposed WENO scheme. Finally, concluding remarks are in Section .

2 The fundamentals of WENO scheme

This part accounts to the construction of the seventh order weighted essentially non-oscillatory scheme in a finite difference framework for approximate the solution of hyperbolic conservation laws

(1)

with initial condition

(2)

Here is a dimensional vector of conserved variables defined for space and time variables respectively, is a flux function which depends on the conserved quantity The system is called hyperbolic if all the eigen values of the Jacobian matrix of the flux function are real and the set of right eigen vectors are complete.

For numerical approximation the spatial domain is discretized with uniform grid, for brevity in the presentation. Let be the length of the cell with center here are known as cell interfaces. The approximation of the spatial derivative in the hyperbolic conservation laws yields a semi-discrete formulation

(3)

Here is an approximation to at a point in time i.e., for the value and are the numerical fluxes which are Lipschitz continuous in each of its arguments and are consistent with the physical flux, that means The conservation property is retrieved by defining a function implicitly through the equation (see Lemma of [32])

(4)

Differentiating (4) with respect to yields

thus should be an approximation to the numerical flux such that

represents the number of cells in a stencil. Thus the numerical flux can be acquired by using higher order polynomial interpolation to with the help of known values of at the cell centers,

To ensure the numerical stability and to avoid entropy violating solutions, the flux is splitted into two parts and the positive and negative parts of respectively, such that

(5)

where and The numerical fluxes and evaluated at reduces (5) as

We will describe here how can be approximated, as is symmetric to the positive part with respect to In the description for the approximation of to follow we drop the sign in the superscript, for simplicity.

2.1 Seventh order WENO scheme

WENO scheme prefers points global stencil, to achieve order of accuracy. The stencil is subdivided into sub-stencils with each sub-stencil bearing cells. In particular, seventh-order WENO scheme accounts to a points stencil, which is subdivided into four -points sub-stencils In accordance with cell each sub-stencil encloses four grid points, specified as

A third degree interpolating polynomial is formulated in each sub-stencil and evaluating it at the cell boundary we retain

(6)

where the coefficients are the Lagrange interpolation coefficients, independent of the values of the flux function but depends on the left shift parameter The equation (6) on each stencil takes the form

The fluxes can be fetched through shifting the index to the left by one in (6). The Taylor’s expansion of the fluxes settle in as

where is the leading order coefficient in the expansion and is the derivative of at The convex combination of these flux functions leads to the approximation of that is, we set

(7)

Here are non-linear weights, satisfying the conditions

If the function is free from discontinuities in all of the sub-stencils we can assess the constants such that the linear combination of provides the seventh order convergence to that is,

The are termed as the ideal weights since they invokes the upstream central scheme of seventh-order, in seven points stencil. The values of these ideal weights are evaluated as

(8)

When the non-linear weights are equal to the ideal weights we have

thus the approximation to the spatial derivative of the flux at the cell-centre is Hence the sufficient condition to achieve the seventh order convergence for the scheme is given by

(9)

where the superscripts and on corresponds to their use in either and respectively. Note that the condition can be weakened for WENO-JS scheme to for a locally Lipschitz continuous function under a suitable condition, which depends on the value of

3 The numerical scheme, WENO-NS7

The essential ingredient of the WENO schemes is in the computation of smoothness indicators, Ha et. al. [11] have established the smoothness indicators measured in norm for the fifth-order WENO scheme. Here we are extending it to the higher order schemes, particularly for the seventh order WENO scheme and calling the scheme as WENO-NS7. The primary thought is of getting a higher order approximation to the derivatives, as it is known that a smoothness indicators based on -norm may give a loss of accuracy in smooth regions [26]. The heuristic construction of smoothness indicators is as follows.

Define the operators which measures the regularity of the solution in each of the four points stencil by estimating the approximate magnitude of derivatives. Once obtained these operators, the smoothness indicators are designed as

(10)

here are the parameters introduced to balance the tradeoff between the accuracy around the smooth region to that of the discontinuous region. The operators are the generalized undivided differences of which are formulated as

(11)

where the coefficient vector

(12)

is obtained by solving the linear system

with . This linear system can be written in the matrix form

with the matrices and defined by

where is the Kronecker delta and there exist a unique solution for this linear system, as is a Vandermode matrix. With the coefficients (12) the operators (11) takes the form

The third operator is the same as in the WENO-BS scheme which is described in appendix. However WENO-NS7 scheme uses the absolute values where as the WENO-BS uses the squared ones. The advantages with these operators is that the approximation of the derivative to be of higher order accuracy at the evaluation point, which is stated in the following theorem.
Thereom Let the stencil and assume that where is an open interval containing . For each the operator (11) satisfies

proof The proof follows in a similar lines that of the Theorem 3.2 of [11]. The Taylor’s expansion of the operators for each and reveals

Thus the operators are in tuned with the Theorem 3.1 stated above.

Now the non-linear weights for the scheme are defined as

(13)

where a global smoothness indicator is taken as

(14)

and the ideal weights, given in (8). To avoid the scenario of zero division a small number is incorporated in the calculations of non-linear weights (13).

Next we discuss the convergence analysis of the WENO-NS7 scheme, in-particularly at the critical points, i.e., we analyze how the weights approaches to the ideal weights in the presence of critical points.

3.1 Convergence order of WENO-NS7 scheme

First consider that there are no critical points and let’s take in (10), from Taylor’s series expansion

(15)

Similarly from the expansion of the global smoothness indicator (14), we’ve

(16)

Then there exist a constant such that

(17)

where .

If , from (15) the smoothness indicators are of the form

then there exists a constant such that

(18)

where . Similarly if then there exist a constant such that

(19)

From (17) and (18), the weights satisfies the sufficient condition (9) if the first derivative vanishes but not the second derivative. In order to satisfy the sufficient condition even at higher order critical points, the nonlinear weights are defined as

(20)

where can be chosen such that the sufficient condition have to hold. Note that from the expansions given in (15-19) and from (20) we have

Clearly, the sufficient condition (9) is satisfies for the WENO-NS7 scheme under the following conditions:

  1. if or if

  2. if .

For numerical verification, value is taken as With these developments, in the next section we’ll test the WENO-NS7 scheme for various examples.

4 Numerical results

Let’s denote the system (1) by

where is an approximation to the derivative In section 2, we have obtained higher order reconstruction for the flux function which is defined in (3). To evolve the solution in time, strong-stability-preserving Runge–Kutta algorithm is used, whose detailed description is in article [10]. The choice of this time integration is to ensure that the order of accuracy for the time evaluation matches with that of the spatial order of accuracy. For linear problems, stage linear method, which is of order is used in the following examples. The method is

For the seventh order the coefficients are given as

This method will not attain for nonlinear problems. So, for nonlinear problems a fourth order nonlinear version of is used with the stability condition where Here is the maximum propagation speed in at time level . The fourth order method is given as

For numerical comparison of the proposed WENO-NS7 scheme, the numerical results are presented along with the numerical results of seventh order WENO-BS [2] and seventh order WENO-Z [5] schemes in the following sections. These results are mostly comparable with WENO-BS scheme in compare to other seventh-order WENO schemes. For completeness these schemes are briefly described in the appendix.

4.1 Scalar test problems

To verify the numerical accuracy and convergence of the scheme examples pertaining to transport and Burger’s equations with various initial profiles are considered. Some of these initial profiles contain jump discontinuity and in some cases, the solution in time leads to shocks. Lax-Friedrich’s flux splitting technique is used in (5). For the scheme WENO-BS, and for WENO-Z and WENO-NS7 schemes is taken along with the CFL number . The parameters in (10) are fixed as and for linear test cases.

4.1.1 Behaviors of new weights

Example 1: For linear advection equation

let the initial condition be

(21)

which is a piecewise continuous function with jump discontinuity at The behavior of the smoothness indicators and the global smoothness indicator for initial time, is displayed in figure 1 for the proposed scheme WENO-NS7.

Figure 1: The values of the smoothness indicators to initial data ( by WENO-NS7 scheme

The approximate solution is computed with uniform discretization of the spatial domain with the step size, up to time with the scheme WENO-NS7 along with WENO-BS and WENO-Z schemes, these solutions are plotted in figure 2 against the exact solution. It can be observed from the plot that the proposed scheme performs better than other schemes near the jump discontinuity. Figure 3 displays the behavior of the weights along with the ideal weights

Figure 2: Numerical solution of linear advection equation with the initial condition (21)
Figure 3: The distribution of ideal weights and the weights .

4.1.2 Accuracy, at critical points

Example 2: Consider the linear transport equation

(22)

with the periodic boundary conditions up to time Three different initial conditions are considered, each of them is a special case to test the convergence analysis.
Case 1: The smooth initial condition

(23)

is taken in this case to verify the order of convergence.

In table 1, the and errors along with the numerical order of accuracy are given for WENO-BS, WENO-Z and WENO-NS7 schemes. It has been observed that the proposed scheme attains the desired order of accuracy and very much efficient than WENO-BS scheme but has the same accuracy as WENO-Z scheme.

WENO-BS WENO-Z WENO-NS7 N 10 3.9281e-03 6.5177e-04 9.9734e-04 20 8.5335e-05 5.52 4.2296e-06 7.26 5.6871e-06 7.45 40 1.2222e-06 6.12 3.3547e-08 6.98 3.3551e-08 7.40 80 1.8418e-08 6.05 2.6304e-10 6.99 2.6304e-10 6.99 160 2.7931e-10 6.04 2.0638e-12 6.99 2.0637e-12 6.99 10 8.5830e-03 1.4372e-03 2.0586e-03 20 2.2095e-04 5.28 6.7946e-06 7.72 1.5028e-05 7.09 40 5.3163e-06 5.38 5.2546e-08 7.01 6.8261e-08 7.78 80 1.4576e-07 5.18 4.1283e-10 6.99 4.3546e-10 7.29 160 4.0863e-09 5.15 3.2415e-12 6.99 3.2736e-12 7.05

Table 1: errors of example (22) with initial condition (23)

Case 2: In this case, the initial condition is chosen as

(24)

which contains first-order critical point, that is, but note that

The and errors along with the numerical order of accuracy are provided in table 2 for WENO-BS, WENO-Z and WENO-NS7 schemes. The proposed scheme WENO-NS7 achieves the desired order of accuracy.

WENO-BS WENO-Z WENO-NS7 N 10 2.9823e-02 —— 3.0006e-02 —– 1.3986e-02 —– 20 7.5922e-04 5.29 5.0945e-04 5.88 3.0947e-04 5.49 40 8.3317e-06 6.51 2.6348e-06 7.59 2.6567e-06 6.86 80 6.9845e-08 6.89 2.1488e-08 6.93 2.1524e-08 6.94 160 7.1684e-10 6.60 1.6933e-10 6.99 1.6934e-10 6.99 10 7.1325e-02 —— 7.0754e-02 —— 3.6883e-02 —— 20 2.2429e-03 4.99 1.3617e-03 5.70 9.2089e-04 5.32 40 3.6731e-05 5.93 7.9578e-06 7.41 7.9697e-06 6.85 80 5.8755e-07 5.96 6.6439e-08 6.90 6.6272e-08 6.91 160 1.1019e-08 5.73 5.2713e-10 6.97 5.2711e-10 6.97

Table 2: errors of example (22) with initial condition (24)

Case 3: The initial condition

(25)

has the nature, but

In table 3, the and errors are tabulated along with the numerical order of accuracy for the schemes WENO-BS, WENO-Z and WENO-NS7. In this case too, the proposed scheme has the desired order of convergence.

WENO-BS WENO-Z WENO-NS7 N 10 2.1618e-01 —— 2.1128e-01 —– 1.9695e-01 —– 20 1.8858e-02 3.52 8.4032e-03 4.65 8.6604e-03 4.50 40 4.1300e-04 5.51 6.4354e-05 7.02 8.6452e-05 6.64 80 5.9000e-06 6.13 4.2948e-07 7.22 4.1937e-07 7.68 160 6.6306e-08 6.47 3.3639e-09 6.99 3.3582e-09 6.96 10 2.9623e-01 —— 2.9814e-01 ——- 3.2319e-01 20 3.8372e-02 2.95 1.4710e-02 4.34 1.7309e-02 4.22 40 1.0319e-03 5.21 1.4691e-04 6.64 2.1456e-04 6.33 80 2.3364e-05 5.46 6.6758e-07 7.78 1.0238e-06 7.71 160 6.1191e-07 5.25 5.2805e-09 6.98 5.5897e-09 7.51

Table 3: errors of example (22) with initial condition (25)

Example 3: Consider the linear advection equation on an interval with the initial condition

(26)

where , , , , , and