The Approximated Semi-Lagrangian WENO Methods Based on Flux Vector Splitting for Hyperbolic Conservation Laws

Fuxing Hu

Department of Mathematics, Huizhou University, Huizhou, China.

DOI: 10.4236/ajcm.2017.71004                

Abstract

The paper is devised to combine the approximated semi-Lagrange weighted essentially non-oscillatory scheme and flux vector splitting. The approximated finite volume semi-Lagrange that is weighted essentially non-oscillatory scheme with Roe flux had been proposed. The methods using Roe speed to construct the flux probably generates entropy-violating solutions. More seriously, the methods maybe perform numerical instability in two-dimensional cases. A robust and simply remedy is to use a global flux splitting to substitute Roe flux. The combination is tested by several numerical examples. In addition, the comparisons of computing time and resolution between the classical weighted essentially non-oscillatory scheme (WENOJS-LF) and the semi-Lagrange weighted essentially non-oscillatory scheme (WENOEL-LF) which is presented (both combining with the flux vector splitting).

Keywords

Semi-Lagrangian Method, WENO Scheme, Flux Splitting

1. Introduction

The semi-Lagrangian methods are popularly used in weather prediction [1] [2] and simulation of the Vlasov equations [3] [4] [5] [6] , and so on. These methods solve the problems with a characteristic tracing algorithm. It is this property of algorithm that leads to the following two advantages: no time discretization and the alleviation of CFL time step restriction. The semi-Lagrangian methods (combined wih WENO reconstruction [7] [8] ) used to solve hyperbolic problems are presented in [9] [10] [11] .

In [10] , the authors proposed a finite volume semi-Lagrangian WENO scheme for advection problems. The scheme combined the Eulerian-Lagrangian framework [12] [13] with high-order WENO reconstruction, and used a integral- based WENO reconstruction to handle trace-back integration. In this framework, the scheme trace each computational Eulerian grid cell at time $t^{n+1}$ backward over a time step along the characteristic line to its Lagrangian track-back region at $t^n$ . The average mass is simply transported from the Lagrangian region to the computational grid cell. In addition, they presented a theoretical proof for the accuracy of the method.

In [11] , the authors developed a finite volume semi-Lagrangian WENO scheme for nonlinear conservation laws. This method can be regarded as the extension of the method for advection problems in [10] . A problem appeared in such extension is that one does not have the particle velocity since it is nonlinearly related to the unknown solution. Hence, one cannot find the exact tracking line of fluid particle. Instead of trying to find the the exact characteristic line of particle, they use a known velocity for the tracing line computation. Since this will not give the correct tracing line, a flux correction procedure is needed to increase the accuracy and numerical stability. The proposed method arrives optimal order of accuracy. However, the procedure of flux correction makes the scheme rather cost to be implemented.

In [9] , a sort of finite difference semi-Lagrangian ENO and WENO schemes are devised for advection problems in incompressible flow. The key part of this paper is that the integral form of advection equation is taken over a triangle region. This integral procedure transforms the flux integration in time at cell center into the integration of mass $u\left(x\mathrm{,}{t}^n\right)$ in space

${\displaystyle {\int }_{t_n}^{t_{n+1}}}f\left(u\left(x_i,t\right)\right)\text{d}t={\displaystyle {\int }_{x_i^{\ast }}^{x_i}}u\left(x,t^n\right)\text{d}x,$

where $x_i^{\ast }$ is backward characteristic point of cell center $x_i$ . This scheme is difficult to extent to nonlinear cases. In addition, there is no optimal linear weights for WENO schemes for solving advection equations with variable coefficients (or also for nonlinear cases).

In [14] , the authors developed an approximated finite volume semi-Lagrangian WENO schemes (WENOEL-Roe) for 1D and 2D nonlinear hyperbolic prolems. The scheme integrates the hyperbolic equation over the control volume

$\left[x_{i-1/2},x_{i+1/2}\right]\times \left[t^n,t^{n+1}\right]$ to obtain a integral equation. They try to directly evaluate the integration of flux function in time at cell edge $x_{i+1/2}$ .

${\displaystyle {\int }_{t_n}^{t_{n+1}}}f\left(u\left(x_{i+\mathrm{1\; <mo>\; /\; </mo>\; 2}}\mathrm{,}t\right)\right)\text{d}t\mathrm{.}$

For linear cases, the integration of flux function in time can be transformed into the integration of interpolation polynomial of flux average in space. For nonlinear cases, a local freeze method is used to freeze the nonlinear into linear cases. The scheme here is different with the traditional semi-Lagrangian scheme [10] [11] . The backward tracing characteristic point is needed in the procedure of evaluating the integration of flux function. The advantages of the WENOEL- Roe scheme are easy implement and high efficiency. Refer to [14] for a detail.

The procedure of evaluating the integration of flux function also depends on the direction of upwind. The Roe average speed is used to identified the upwinding. It is known that the Roe schemes maybe generates entropy-violating solutions. More seriously, these scheme can perform numerical instability for some two-dimensional problems [15] [16] [17] . A local entropy correction can be used to remedy this deficiency. However, it is usually more robust to use a global flux splitting. In this paper, an approximated finite volume semi- Lagrangian WENO scheme with the smooth Lax-Friedrichs flux splitting (WENOEL-LF) is presented. The WENOEL-LF scheme is less resolution than WENOEL-Roe since the Lax-Friedrichs flux splitting method is more dissipative than Roe method. However, the advantage of WENOEL-LF scheme is it generally generates the entropy solution.

In the paper that follows, we will review the WENOEL-Roe scheme briefly in Section 2. In Section 3, we will present the formulation of WENOEL-LF scheme in Section 3 in detail. The comparisons of resolution and computing time between the WENOEL-LF and WENOJS-LF schemes is presented in Section 4.

2. Review the WENOEL-Roe Scheme

Consider the linear advection equation

$u_t\left(x,t\right)+f_x\left(u\left(x,t\right)\right)=0,a\le x\le b,t\ge 0,$ (1)

where $f_x\left(u\left(x,t\right)\right)=au\left(x,t\right)$ . Integrating the Equation (1) on control volume $\left[x_{i-1/2}\mathrm{,}{x}_{i+1/2}\right]\times \left[t^n\mathrm{,}{t}^{n+1}\right]$ and rearranging this equation gives

$\begin{array}{c}\frac{1}{\Delta x}{\displaystyle {\int }_{x_{i-1/2}}^{x_{i+1/2}}}u\left(x,t_{n+1}\right)\text{d}x=\frac{1}{\Delta x}{\displaystyle {\int }_{x_{i-1/2}}^{x_{i+1/2}}}u\left(x,t_n\right)\text{d}x+\frac{\Delta t}{\Delta x}(\frac{1}{\Delta t}{\displaystyle {\int }_{t^n}^{t^{n+1}}}f\left(u\left(x_{i-1/2},t\right)\right)\text{d}t\\ -\frac{1}{\Delta t}{\displaystyle {\int }_{t^n}^{t^{n+1}}}f\left(u\left(x_{i+1/2},t\right)\right)\text{d}t).\end{array}$ (2)

Denote the $i$ -th cell average by

${\bar{u}}_i^n=\frac{1}{\Delta x}{\displaystyle {\int }_{x_{i-1/2}}^{x_{i+1/2}}}u\left(x,t_n\right)\text{d}x,$

and average flux at cell edge $x_{i+1/2}$ by

${\bar{f}}_{i+1/2}^n=\frac{1}{\Delta t}{\displaystyle {\int }_{t_n}^{t_{n+1}}}f\left(u\left(x_{i+1/2},t\right)\right)\text{d}t,$ (3)

then (2) can be written as

${\bar{u}}_i^{n+1}={\bar{u}}_i^n-\frac{\Delta t}{\Delta x}\left({\bar{f}}_{i+1/2}^n-{\bar{f}}_{i-1/2}^n\right).$ (4)

Denote the value $U_i^n$ approximates the average value ${\bar{u}}_i^j$ and the numerical flux $F_{i+1/2}^n$ approximates flux function ${\bar{f}}_{i+1/2}^n$ . Then we obtain

$U_i^{n+1}=U_i^n-\frac{\Delta t}{\Delta x}\left(F_{i+1/2}^n-F_{i-1/2}^n\right).$ (5)

For evaluating average flux ${\bar{f}}_{i+1/2}^n$ in (3), we can firstly apply a 5th-order reconstruction based on piecewise constant average fluxes ${\left\{f_k\right\}}_{k=1}^N$ on stencil $S_i=\left\{I_{i-2},I_{i-1},I_i,I_{i+1},I_{i+2}\right\}$ to obtain interpolation polynomial $P_i\left(x\right)$ on cell $I_i$

$f\left(u\left(x,t_n\right)\right)=P_i\left(x\right)+\mathcal{O}\left(\Delta x^5\right)\mathrm{.}$ (6)

Here, we assume the advection velocity $a>0$ , then $u\left(x\mathrm{,}{t}_n\right)$ simply spreads to the right with velocity $a$ , which gives

$f\left(u\left(x_{i+1/2},t\right)\right)=f\left(u\left(x_{i+1/2}-a\left(t-t_n\right),t_n\right)\right)$ (7)

in any time when $t\in \left[t_n\mathrm{,}{t}_{n+1}\right]$ . Combining the formulas (6) with (7), the flow rate $f\left(u\left(x\mathrm{,}t\right)\right)$ at cell edge $x_{i+1/2}$ is

$f\left(u\left(x_{i+1/2},t\right)\right)=P_i\left(x_{i+1/2}-a\left(t-t_n\right)\right)+\mathcal{O}\left(\Delta x^5\right)\mathrm{.}$ (8)

In this case, the average flux ${\bar{f}}_{i+1/2}^n$ can be expressed approximately as

$\begin{array}{c}{\bar{f}}_{i+1/2}^n=\frac{1}{\Delta t}{\displaystyle {\int }_{t_n}^{t_{n+1}}}f\left(u\left(x_{i+1/2}\mathrm{,}t\right)\right)\text{d}t\\ =\frac{1}{\Delta t}{\displaystyle {\int }_{t_n}^{t_{n+1}}}{P}_i\left(x_{i+1/2}-a\left(t-t_n\right)\right)\text{d}t+\mathcal{O}\left(\Delta x^5\right)\\ =\frac{1}{a\Delta t}{\displaystyle {\int }_{x_{i+1/2}-a\Delta t}^{x_{i+1/2}}}{P}_i\left(x\right)\text{d}x+\mathcal{O}\left(\Delta x^5\right)\mathrm{.}\end{array}$ (9)

Omitting the high-order term $\mathcal{O}\left(\Delta x^5\right)$ , we obtain the numerical flux

$F_{i+1/2}^n\mathrm{<mo>\; =\; </mo>}\frac{1}{a\Delta t}{\displaystyle {\int }_{x_{i+1/2}-a\Delta t}^{x_{i+1/2}}}{P}_i\left(x\right)\text{d}x.$ (10)

The last equality in (9) is obtained by integration of substitution

$x=x_{i+1/2}-a\left(t-t_n\right)$ . From the Equation (9), the integration in time $\left[t_n\mathrm{,}{t}_{n+1}\right]$ is transformed into integration in space $\left[x_{i+1/2}-a\Delta t\mathrm{,}{x}_{i+1/2}\right]$ . Due to the integrand $P_i\left(x\right)$ is reconstruction polynomial, the last integration in (10) can be computed exactly.

In the following of this section, we will present the 5th-order WENO reconstruction process to approximate the integral in (10)

${\displaystyle {\int }_{x_{i+1/2}-a\Delta t}^{x_{i+1/2}}}{P}_i\left(x\right)\text{d}x\mathrm{.}$

The 5th-order WENO reconstruction procedure is represented as the convex combination of three 3rd-order reconstructions. First, we intend to reconstruct three 3rd-order conservative polynomials $P_i^j\left(x\right)$ on cell $I_i=\left[x_{i-1/2},x_{i+1/2}\right]$ based on the piecewise constant average fluxes on stencils

$S_i^0=\left\{I_i,I_{i+1},I_{i+2}\right\},S_i^1=\left\{I_{i-1},I_i,I_{i+1}\right\},S_i^2=\left\{I_{i-2},I_{i-1},I_i\right\},$

where the subscript $i$ denotes the polynomial on cell $I_i$ and superscript $j$ denotes the reconstruction based on stencil $S_i^j$ . So far, we have obtained the conservative interpolation polynomials $P_i^j\left(x\right)$ on each stencil $S_i^j$ . In the end, the integral in (10) can be expressed as

${\displaystyle {\int }_{x_{i+1/2}-a\Delta t}^{x_{i+1/2}}}{P}_i\left(x\right)\text{d}x={\displaystyle \underset{j\mathrm{<mo>\; =\; </mo>\; 0}}{\overset{2}{\sum }}}{d}^j{\displaystyle {\int }_{x_{i+1/2}-a\Delta t}^{x_{i+1/2}}}{P}_i^j\left(x\right)\text{d}x\mathrm{,}$

where $d^j$ are the optimal weights. To alleviate the effect of the non-smooth stencils, the nonlinear weights can be constructed as follows

$w_i^j=\frac{a_i^j}{a_i^0+a_i^1+a_i^2}$ (11)

and

$a_i^j=\frac{d^j}{{\left(\epsilon +{\beta }_i^j\right)}^2},$ (12)

where ${\beta }_i^j$ is the indicator of smoothness of the polynomial on the stencil $S_i^j$ . Finally, the numerical flux should be expressed as

$F_{i+1/2}^n=\frac{1}{a\Delta t}{\displaystyle \underset{j=0}{\overset{2}{\sum }}}{w}_i^j{\displaystyle {\int }_{x_{i+1/2}-a\Delta t}^{x_{i+1/2}}}{P}_i^j\left(x\right)\text{d}x\mathrm{.}$ (13)

Substituting the formula of polynomial $P_i^j\left(x\right)$ into (13) gives

$\begin{array}{c}{F}_{i+1/2}^n=\frac{1}{a\Delta t}{\displaystyle \underset{j=0}{\overset{2}{\sum }}}{w}_i^j{\displaystyle {\int }_{x_{i+1/2}-a\Delta t}^{x_{i+1/2}}}{P}_i^j\left(x\right)\text{d}x\\ =\frac{1}{6}[w_i^0\left({\lambda }^2\left(f_i-2f_{i+1}+f_{i+2}\right)+\lambda \left(3f_i-3f_{i+1}\right)+\left(2f_i+5f_{i+1}-f_{i+2}\right)\right)\\ +w_i^1\left({\lambda }^2\left(f_{i-1}-2f_i+f_{i+1}\right)+\lambda \left(3f_i-3f_{i+1}\right)+\left(-f_{i-1}+5f_i+2f_{i+1}\right)\right)\\ +w_i^2\left({\lambda }^2\left(f_{i-2}-2f_{i-1}+f_i\right)+\lambda \left(-3f_{i-2}+9f_{i-1}-6f_i\right)+\left(2f_{i-2}-7f_{i-1}+11f_i\right)\right)].\end{array}$ (14)

In contrast, when the advection velocity $a<0$ , $u\left(x\mathrm{,}{t}_n\right)$ spreads to the left with velocity $a$ , which gives

$f\left(u\left(x_{i+1/2},t\right)\right)=f\left(u\left(x_{i+1/2}-a\left(t-t_n\right),t_n\right)\right)$ (15)

and the flow rate $f\left(u\left(x\mathrm{,}t\right)\right)$ at cell edge $x_{i+1/2}$ is

$f\left(u\left(x_{i+1/2}\mathrm{,}t\right)\right)=P_{i+1}\left(x_{i+\mathrm{1\; <mo>\; /\; </mo>\; 2}}-a\left(t-t_n\right)\right)+\mathcal{O}\left(\Delta x^5\right)\mathrm{.}$ (16)

Similarly, the average flux $F_{i+1/2}^n$ can be expressed as

$F_{i+1/2}^n=\frac{1}{-a\Delta t}{\displaystyle \underset{j=0}{\overset{2}{\sum }}}{w}_{i+1}^j{\displaystyle {\int }_{x_{i+1/2}}^{x_{i+1/2}-a\Delta t}}{P}_{i+1}^j\left(x\right)\text{d}x.$ (17)

The nonlinear weight $w_{i+1}^j$ can be constructed similarly as (11) and (12). Substituting the formula of polynomial $P_{i+1}^j\left(x\right)$ into (17) gives

$\begin{array}{c}{F}_{i+1/2}^n=\frac{1}{6}[w_{i+1}^0({\lambda }^2\left(f_{i+1}-2f_{i+2}+f_{i+3}\right)+\lambda \left(6f_{i+1}-9f_{i+2}+3f_{i+3}\right)\\ +\left(11f_{i+1}-7f_{i+2}+2f_{i+3}\right))+w_{i+1}^1({\lambda }^2\left(f_i-2f_{i+1}+f_{i+2}\right)\\ +\lambda \left(3f_i-3f_{i+1}\right)+\left(2f_i+5f_{i+1}-f_{i+2}\right))\\ +w_{i+1}^2\left({\lambda }^2\left(f_{i-1}-2f_i+f_{i+1}\right)+\lambda \left(3f_i-3f_{i+1}\right)+\left(-f_{i-1}+5f_i+2f_{i+1}\right)\right)]\mathrm{.}\end{array}$ (18)

The above process of approximating average flux is reasonable for linear advection equation. For nonlinear problems, the formulas (9) and (17) is not hold any more since the solution no longer simply translates uniformly. And generally the tracking back points cannot be found exactly (even cannot find the points with high accuracy). Hence for nonlinear case, rather than trying to find the tracking back points; we freeze the nonlinear equation to linear formation locally and apply the procedure above to it. For solving the nonlinear case, the propagation direction is distinguished by Rankine-Hugoniot jump conditions

and propagation velocity is chosen to be $a={\frac{\partial f\left(u\right)}{\partial u}|}_{u=u_i}$

( $\text{if}\text{propagation}\text{direction}>0$ ) or $a={\frac{\partial f\left(u\right)}{\partial u}|}_{u=u_{i+1}}$

( $\text{if}\text{propagation}\text{direction}<0$ ).

3. Formulation of WENOEL-LF Scheme

In this section, we solve nonlinear problems by using a more robust global flux splitting

$f\left(u\right)=f^{+}\left(u\right)+f^{-}\left(u\right)$ (19)

where

$\frac{\text{d}{f}^{+}\left(u\right)}{\text{d}u}\ge 0\text{and}\frac{\text{d}{f}^{-}\left(u\right)}{\text{d}u}\le 0.$

Insert the Equation (19) into (1) and (2), we can obtain

${\bar{u}}_i^{n+1}={\bar{u}}_i^n-\frac{\Delta t}{\Delta x}\left({\bar{f}}_{i+1/2}^{n,+}-{\bar{f}}_{i-1/2}^{n,+}\right)-\frac{\Delta t}{\Delta x}\left({\bar{f}}_{i+1/2}^{n,-}-{\bar{f}}_{i-1/2}^{n,-}\right),$ (20)

which is similar to conservation formula (4), and denote

${\bar{f}}_{i+1/2}^{n,\pm }=\frac{1}{\Delta t}{\displaystyle {\int }_{t^n}^{t^{n+1}}}{f}^{\pm }\left(u\left(x_{i+1/2},t\right)\right)\text{d}t\mathrm{.}$

A scheme approximated (20) can be written as

$U_i^{n+1}=U_i^n-\frac{\Delta t}{\Delta x}\left(F_{i+1/2}^{n,+}-F_{i-1/2}^{n,+}\right)-\frac{\Delta t}{\Delta x}\left(F_{i+1/2}^{n,-}-F_{i-1/2}^{n,-}\right).$ (21)

This scheme is conservative. Since if we sum $U_i^{n+1}$ over the whole set of cells, we obtain

${\displaystyle \underset{i=1}{\overset{N}{\sum }}}{U}_i^{n+1}={\displaystyle \underset{i=1}{\overset{N}{\sum }}}{U}_i^n-\frac{\Delta t}{\Delta x}\left(\left(F_{N+1/2}^{n,+}+F_{N+1/2}^{n,-}\right)-\left(F_{1/2}^{n,+}+F_{1/2}^{n,-}\right)\right)\mathrm{.}$

The formulas $F_{N+1/2}^{n,+}+F_{N+1/2}^{n,-}$ and $F_{1/2}^{n,+}+F_{1/2}^{n,-}$ denote the fluxes at the extreme edges. The sum of the flux differences cancels out except for the fluxes at the extreme.

The simplest smooth flux splitting we chose is the Lax-Friedrichs splitting

$f^{+}\left(u\right)=\frac{1}{2}\left(f\left(u\right)+\alpha u\right),f^{-}(u)=\frac{1}{2}\left(f\left(u\right)-\alpha u\right),$

where $\alpha$ is taken as $\alpha =\underset{1\le u_i\le N}{{\displaystyle \max }}\left|\frac{\text{d}f\left(u\right)}{\text{d}u}\right|$ over the whole set of cell averages. Since $\frac{\text{d}{f}^{+}\left(u\right)}{\text{d}u}\ge 0$ , the fluxes $f^{+}\left(u\right)$ flow right through the cell edges. The flow

velocity $a^{+}$ for flux $f^{+}\left(u\right)$ at cell edge $x_{i+1/2}$ is chosen as

$a^{+}=\frac{1}{2}\left({\frac{\text{d}{f}^{+}\left(u\right)}{\text{d}u}|}_{u=u_i}+\alpha \right)$ . Similarly, The flow velocity $a^{-}$ for flux $f^{-}\left(u\right)$ at cell edge $x_{i+1/2}$ is chosen as $a^{-}=\frac{1}{2}\left({\frac{\text{d}{f}^{-}\left(u\right)}{\text{d}u}|}_{u=u_{i+1}}-\alpha \right)$ .

Inserting $f^{+}\left(u\right)$ and $f^{-}\left(u\right)$ into (14) and (18), respectively, the numerical fluxes can be obtained as follows,

$\begin{array}{c}{F}_{i+1/2}^{n\mathrm{,}+}=\frac{1}{6}[w_i^0\left({\left({\lambda }^{+}\right)}^2\left(f_i^{+}-2f_{i+1}^{+}+f_{i+2}^{+}\right)+{\lambda }^{+}\left(3f_i^{+}-3f_{i+1}^{+}\right)+\left(2f_i^{+}+5f_{i+1}^{+}-f_{i+2}^{+}\right)\right)\\ +w_i^1({\left({\lambda }^{+}\right)}^2\left(f_{i-1}^{+}-2f_i^{+}+f_{i+1}^{+}\right)+{\lambda }^{+}\left(3f_i^{+}-3f_{i+1}^{+}\right)\\ +\left(-f_{i-1}^{+}+5f_i^{+}+2f_{i+1}^{+}\right))+w_i^2({\left({\lambda }^{+}\right)}^2\left(f_{i-2}^{+}-2f_{i-1}^{+}+f_i^{+}\right)\\ +{\lambda }^{+}\left(-3f_{i-2}^{+}+9f_{i-1}^{+}-6f_i^{+}\right)+\left(2f_{i-2}^{+}-7f_{i-1}^{+}+11f_i^{+}\right))],\end{array}$ (22)

where ${\lambda }^{+}=\frac{a^{+}\Delta t}{\Delta x}$ , and

$\begin{array}{c}{F}_{i+1/2}^{n,-}=\frac{1}{6}[w_{i+1}^0({\left({\lambda }^{-}\right)}^2\left(f_{i+1}^{-}-2f_{i+2}^{-}+f_{i+3}^{-}\right)+{\lambda }^{-}\left(6f_{i+1}^{-}-9f_{i+2}^{-}+3f_{i+3}^{-}\right)\\ +\left(11f_{i+1}^{-}-7f_{i+2}^{-}+2f_{i+3}^{-}\right))+w_{i+1}^1({\left({\lambda }^{-}\right)}^2\left(f_i^{-}-2f_{i+1}^{-}+f_{i+2}^{-}\right)\\ +{\lambda }^{-}\left(3f_i^{-}-3f_{i+1}^{-}\right)+\left(2f_i^{-}+5f_{i+1}^{-}-f_{i+2}^{-}\right))\\ +w_{i+1}^2\left({\left({\lambda }^{-}\right)}^2\left(f_{i-1}^{-}-2f_i^{-}+f_{i+1}^{-}\right)+{\lambda }^{-}\left(3f_i^{-}-3f_{i+1}^{-}\right)+\left(-f_{i-1}^{-}+5f_i^{-}+2f_{i+1}^{-}\right)\right)]\end{array}$ (23)

where ${\lambda }^{-}=\frac{a^{-}\Delta t}{\Delta x}$ .

Algorithm

Here, we conclude the algorithm for computing the approximated solution $U_i^{n+1}$ at time step $t^{n+1}$ .

1) Split the flux function $f\left(u\right)$ into the positive flux $f^{+}\left(u\right)$ and negative flux $f^{-}\left(u\right)$ .

2) Determine the flow velocities $a^{+}$ and $a^{-}$ of the positive flux $f^{+}\left(u\right)$ and negative flux $f^{-}\left(u\right)$ , respectively.

3) Compute the numerical fluxes $F_{i+1/2}^{n,\pm }$ and $F_{i-1/2}^{n\mathrm{,}\pm }$ by formulas (22) and (23).

4) Insert the numerical fluxes computed above into (21), we obtain the approximated solution $U_i^{n+1}$ at time step $t^{n+1}$ .

4. Numerical Results

In this section, we use several 1D and 2D nonlinear examples to test the WENOEL-LF scheme. The comparisons of resolution and computing time between the WENOEL-LF and WENOJS-LF is presented. It is found that, with the same number of cells, the WENOJS-LF scheme has slightly higher resolution than WENOEL-LF scheme. However, the computing time of WENOJS-LF scheme is almost two times for scalar equation (and almost three times for nonlinear system) over that of WENOEL-LF scheme.

4.1. Burgers’ Equation

Consider the inviscid Burgers’ equation

$u_t+{\left(\frac{u^2}{2}\right)}_x=0$ (24)

with two initial conditions:

$u\left(x,0\right)=\{\begin{array}{cc}1& -1\le x\le 0,\\ 0& 0\le x\le 1,\end{array}$ (25)

and

$u\left(x,0\right)=1.5-\sin \left(\text{π}x\right).$ (26)

The Burger’s Equation (24) with discontinuous initial condition (25) develops the solution which consists of a rarefaction wave and a shock wave. The numerical solutions computed by the WENOEL-LF and WENOJS-LF schemes are shown in Figure 1. These two solutions are both computed with $N=100$ and CFL = 0.1. The final output time is chosen to be $t=1$ . From Figure 1, we can find that the solution of WENOJS-LF scheme has slightly better resolution than that of WENOEL-LF scheme, especially around the rarefaction wave. And, in solving nonlinear cases (including the following tests), the WENOJS-LF scheme generally possesses slightly higher resolution than WENOEL-LF scheme when the same amount of cells is used.

Although the WENOJS-LF scheme has higher resolution than WENOEL-LF scheme, the WENOEL-LF scheme has the advantage of decreasing the computing time. Table 1 presents the $l_1$ error and computing time. With the same number of cells, the errors generated by WENOEL-LF are larger than WENOJS-LF scheme, but the computing times of WENOEL-LF are only half of

(a) (b)

Figure 1. The solutions of problem (24) (25) are computed by the WENOEL-LF and WENOJS-LF methods: (a) WENOEL-LF scheme; (b) WENOJS-LF scheme.

that of the WENOJS-LF scheme. For comparing the efficiency between these two schemes, we plots the relationship between the $l_1$ error and computing time for these two schemes. From Figure 2, we can find that the efficiency of WENOEL-LF is higher than that of WENOJS-LF scheme. That is, to achieve the same $l_1$ error, the WENOEL-LF needs less computing time.

For the problem (24) (26), the initial condition is smooth, and the solution

evolves discontinuity at $t=\frac{1}{\text{π}}$ . The Figure 3 is plotted with $N=100$ and output

time $t=0.3$ . From this figure, no obviously difference is presented (however, by closed inspection, the resolution of WENOJS-LF is still slightly better than WENOEL-LF scheme). Similar to the last test, the superiority of our scheme lies in decreasing the computing time. Therefore, in efficiency, the WENOEL-LF scheme is still advantageous over the WENOJS-LF scheme.

Table 1. The comparisons of computing time (in seconds) and $l_1$ error for initial value problem (24) (25) between the WENOEL-LF and WENOJS-LF methods.

Figure 2. The comparison of efficiency for problem (24) (25) is presented between the WENOEL-LF and WENOJS-LF methods.

(a) (b)

Figure 3. The solutions of problem (24) (26) are computed by the WENOEL-LF and WENOJS-LF methods: (a) WENOEL-LF scheme; (b) WENOJS-LF scheme.

4.2. The 1D Euler Equation

In this subsection, we consider 1D Euler equations since one of the main application areas of high-resolution scheme is compressible gas dynamics,

${\left(\begin{array}{c}\rho \\ \rho u\\ E\end{array}\right)}_t+{\left(\begin{array}{c}\rho u\\ \rho u^2+p\\ \left(E+p\right)u\end{array}\right)}_x=0,$ (27)

where $\rho$ , $u$ , $p$ , $E$ are density, velocity, pressure and total energy, respectively. The system of equations is closed by the equation of state for an ideal polytropic gas:

$E=\frac{p}{\gamma -1}+\frac{1}{2}\rho u^2,$

where the ratio of specific heats $\gamma =1.4$ .

The following three initial conditions combined with Euler Equation (27) are considered, which are often used to examine the methods for solving Euler equations:

$\left(\rho ,u,p\right)=\{\begin{array}{ll}\left(0.445,0.698,3.528\right),\hfill & \text{if}-5\le x<0,\hfill \\ \left(0.5,0,0.571\right),\hfill & \mathrm{if}0\le x\le 5,\hfill \end{array}$ (28)

$\left(\rho ,v,p\right)=\{\begin{array}{ll}\left(\frac{27}{7},\frac{4\sqrt{35}}{9},\frac{31}{3}\right),\hfill & \text{if}-5\le x<-4,\hfill \\ \left(1+0.2\sin \left(5x\right),0,1\right),\hfill & \text{if}-4\le x\le 5,\hfill \end{array}$ (29)

$\left(\rho ,u,p\right)=\{\begin{array}{ll}\left(1,0,1000\right),\hfill & \text{if}-5\le x<-4,\hfill \\ \left(1,0,0.01\right),\hfill & \text{if}-4\le x<4,\hfill \\ \left(1,0,100\right),\hfill & \text{if}4\le x\le 5.\hfill \end{array}$ (30)

For the problem (27) (28), called Lax problem [18] , we solve it with

$\Delta t=0.1\Delta x/\lambda$ for the WENOEL-LF and WENOJS-LF methods. The Figure 4 shows the numerical solutions computed by the WENOEL-LF and WENOJS-LF schemes with $N=200$ . From this figure, no obviously difference can be found between these two schemes and they both present high-resolution results. To make further comparison, we show the $l_1$ error and computing time in Table 2. The exact solution of Riemann problem (27) (28) is generated by the code of E. F. Toro in [19] . From this table, we can find that, with the same number of cells, the $l_1$ error of WENOEL-LF is slightly larger than WENOJS-LF scheme. But the WENOEL-LF scheme only need almost one-third of the computing time compared with the WENOJS-LF scheme. Here, we also plot the relationship between the $l_1$ error and computing time in Figure 5. Apparently, to achieve the equal error, the WENOEL-LF scheme spends less computing time than WENOJS-LF scheme.

(a) (b)

Figure 4. The solutions of problem (27) (28) are computed by the WENOEL-LF and WENOJS-LF methods; (a) WENOEL-LF scheme; (b) WENOJS-LF scheme.

Table 2. The comparisons of computing time (in seconds) and $l_1$ error for the problem (27) (28) between the WENOEL-LF and WENOJS-LF methods.

For solving the problem (27) (29), we have the same set as the previous example but with different number of cells. Because no exact solution can be obtained, there is no comparison of error. In Figure 6, we show the numerical

Figure 5. The comparison of efficiency for problem (27) (28) is presented between the WENOEL-LF and WENOJS-LF methods.

Figure 6. The numerical solutions for problem (27) (29) is computed by the WENOEL-LF and WENOJS-LF methods. The “+” and “ $\circ$ ” denotes the solutions computed by the WENOEL-LF and WENOJS-LF schemes with $N=400$ , respectively. The “×” denotes the solutions computed by the WENOEL-LF scheme with $N=600$ .

solution of the WENOJS-LF scheme with $N=400$ and the numerical solutions of the WENOEL-LF scheme with $N=400$ and $N=600$ , respectively. The Figure 7 is the zoomed version of the Figure 6 around the high-frequency wave, From this zoomed figure, it is easily found that, with the same number of cells $N=400$ , the WENOJS-LF scheme generates the solution with higher resolution than WENOEL-LF scheme. However, the WENOJS-LF and WENOEL-LF schemes spend about 8 and 3 seconds to evolve the solutions, respectively. In addition, we present another numerical solution generated by the WENOEL-LF scheme with $N=600$ . This solution spends about 7 seconds and performs clearly higher resolution. In a word, with the same computing time, the WENOEL-LF scheme can obtains the solution with higher resolution.

For the problem (27) (30), also due to there is no exact solution, we use the same way as last example to compare these two schemes. In Figure 8, we show the numerical solution by the WENOJS-LF scheme with $N=400$ and the numerical solutions by the WENOEL-LF scheme with $N=400$ and $N=600$ , respectively. With the same number of cells, as usual, the WENOJS-LF scheme gives higher-resolution solution than the WENOEL-LF scheme, which is shown in Figure 9. The computing times needed for the WENOJS-LF and WENOEL-LF schemes, when $N=400$ , are about 15 and 6 seconds, respectively. If the number of cells increases from $N=400$ to $N=600$ for the WENOEL-LF scheme, the numerical solution has higher resolution than the WENOJS-LF scheme which is computed with $N=400$ . However, it is noted that the WENOEL-LF scheme only needs 13 seconds at this time. That is, to achieve the same resolution the WENOEL-LF scheme needs much more cells, but takes less computing time.

Figure 7. The zoomed version of Figure 6 around the part of high-frequency wave.

Figure 8. The numerical solutions for problem (27) (30) is computed by the WENOEL-LF and WENOJS-LF methods. The “+” and “ $\circ$ ” denotes the solutions computed by the WENOEL-LF and WENOJS-LF schemes with $N=400$ , respectively. The “×” denotes the solutions computed by the WENOEL-LF scheme with $N=600$ .

Figure 9. The zoomed version of Figure 8 around complicated region.

4.3. The 2D Euler System

Finally, we consider a numerical experiment for 2D Euler equations for gas dynamics,

${\left[\begin{array}{c}\rho \\ \rho u\\ \rho v\\ E\end{array}\right]}_t+{\left[\begin{array}{c}\rho u\\ \rho u^2+p\\ \rho uv\\ \left(E+p\right)u\end{array}\right]}_x+{\left[\begin{array}{c}\rho v\\ \rho uv\\ \rho v^2+p\\ \left(E+p\right)v\end{array}\right]}_y=0,$ (31)

where the equation of state is

$p=\left(\gamma -1\right)\left(E-\frac{1}{2}\rho \left(u^2+v^2\right)\right),\gamma =\mathrm{1.4.}$

Here, we apply the WENOEL-LF and WENOJS-LF schemes to the 2D double-Mach shock reflection problem where a strong vertical shock moves horizontally into a wedge which inclined with some angle with the ratio of specific heats $\gamma =1.4$ . Initially, this problem was proposed by Woodward and Colella [20] and had been taken extensively as a test example for high-order schemes. The computational domain is chosen to be $\left[0,4\right]\times \left[0,1\right]$ and the

reflective wall lies on the bottom of the computational domain for $\frac{1}{6}\le x\le 4$ . In the beginning, a Mach 10 shock, moving right, is located at $x=\frac{1}{6}$ , $y=0$ and

makes an angle 60˚ with the x-axis. For the boundary conditions, the exact

postshock condition is imposed for bottom boundary from $x=0$ to $x=\frac{1}{6}$ ,

and the reflective boundary condition is imposed for the rest; the flows are imposed on the top boundary such that there is no interaction with the Mach 10 shock; inflow and outflow boundary conditions are set for the left and right boundaries respectively. The unshocked fluid has a density of 1.4, a pressure of 1 and this problem is run until $t=0.2$ .

In Figure 10, we plot the numerical solution of the WENOJS-LF scheme with cells 400 × 100 and the numerical solutions of the WENOEL-LF scheme with cells 400 × 100 and 560 × 140, respectively. As the conclusion of 1D cases, with the same number of cells, the WENOJS-LF scheme can obtain higher-resolution solutions than the WENOEL-LF scheme, but with much more computing time. When we increase the number of cells from 400 × 100 to 560 × 140, the numerical solution of the WENOEL-LF scheme performs higher resolution but with almost the same computing time as the WENOJS-LF scheme with 400 × 100.

5. Conclusion

This paper aims to combine the finite volume semi-Lagrangian WENO method with flux vector splitting. The proposed scheme is rather robust and easily implemented. And the computing time of the WENOEL-LF scheme is about half and one third of the WENOJS-LF scheme in scalar equation and nonlinear system, respectively. But with the same number of cells, unlike the performance in [14] , the numerical solution of WENOJS-LF scheme is slightly better than that of WENOEL-LF scheme. If we compare the $L_1$ error with computing time (or the resolution of numerical solution and computing time), then we can find the

Figure 10. Double Mach problem. (a) The density contour is computed by the WENOEL-LF scheme with 400 × 100; (b) The density contour is computed by the WENOEL-LF scheme with 560 × 140; (c) The density contour is computed by the WENOJS-LF scheme with 400 × 100. 30 equally spaced contour lines are plotted from 1.731 to 22.9705.

WENOEL-LF scheme possesses higher efficiency than the WENOJS-LF scheme. That is, to obtain the numerical solution with the same resolution, the WENOEL-LF scheme needs less computing time.

Acknowledgements

This work was supported by the National Natural Science Foundation of China (No. 11501238), Natural Science Foundation of Guangdong Province (No. 2012A030313119) and supported by the Major Project Foundation of Guangdong Province Education Department (No. 2014KZDXM070).

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Staniforth, A. and Cote, J. (1991) Semi-Lagrangian Integration Schemes for Atmospheric Models—A Review. Monthly Weather Review, 119, 2206-2223. https://doi.org/10.1175/1520-0493(1991)119<2206:SLISFA>2.0.CO;2
[2]	Lin, S.J. and Rood, R.B. (1996) Multi-Dimensional Flux-Form Semi-Lagrangian Transport Schemes. Monthly Weather Review, 124, 2046-2070. https://doi.org/10.1175/1520-0493(1996)124<2046:MFFSLT>2.0.CO;2
[3]	Sonnendrucker, E., Roche, J., Bertrand, P. and Ghizzo, A. (1999) The Semi-Lagran gian Method for the Numerical Resolution of the Vlasov Equation. Journal of Computational Physics, 149, 201-220.
[4]	Crouseilles, N., Mehrenberger, M. and Sonnendrucker, E. (2010) Conservative Semi-Lagrangian Schemes for Vlasov Equations. Journal of Computational Physics, 229, 1927-1953.
[5]	Qiu, J.M. and Shu, C.W. (2011) Conservative Semi-Lagrangian Finite Difference WENO Formulations with Applications to the Vlasov Equation. Journal of Computational Physics, 10, 979-1000. https://doi.org/10.4208/cicp.180210.251110a
[6]	Qiu, J.M. and Christlieb, A. (2010) A Conservative High Order Semi-Lagrangian WENO Method for the Vlasov Equation. Journal of Computational Physics, 229, 1130-1149.
[7]	Liu, X.D., Osher, S. and Chan, T. (1994) Weighted Essentially Non-Oscillatory Schemes. Journal of Computational Physics, 115, 200-212.
[8]	Jiang, G.S. and Shu, C.W. (1996) Efficient Implementation of Weighted ENO Schemes. Journal of Computational Physics, 126, 202-228.
[9]	Arbogast, T. and Wheeler, M.F. (1995) A Characteristics-Mixed Finite Element Method for Advection-Dominated Transport Problems. SIAM Journal on Numerical Analysis, 32, 404-424. https://doi.org/10.1137/0732017
[10]	Qiu, J.M. and Shu, C.W. (2011) Conservative High Order Semi-Lagrangian Finite Difference WENO Methods for Advection in Incompressible Flow. Journal of Computational Physics, 230, 863-889.
[11]	Arbogast, T. and Huang, C. (2006) A Fully Mass and Volume Conserving Implementation of a Characteristic Method for Transport Problems. SIAM Journal on Scientific Computing, 28, 2001-2022. https://doi.org/10.1137/040621077
[12]	Huang, C., Arbogast, T. and Qiu, J. (2012) An Eulerian-Lagrangian WENO Finite Volume Scheme for Advection Problems. Journal of Computational Physics, 231, 4028-4052.
[13]	Huang, C. and Arbogast, T. (2017) An Eulerian–Lagrangian Weighted Essentially Nonoscillatory scheme for Nonlinear Conservation Laws. Numerical Methods for Partial Differential Equations, 33, 651-680.
[14]	Sanders, R., Morano, E. and Druguet, M. (1998) Multi-Dimensional Dissipation for Upwind Schemes: Stability and Applications to Gas Dynamics. Journal of Computational Physics, 145, 511-537.
[15]	Pandolfi, M. and D’Ambrosio, D. (2001) Numerical Instabilities in Upwind Methods: Analysis and Cures for the “Carbuncle” Phenomenon. Journal of Computational Physics, 166, 271-301.
[16]	Kim, S., Kim, C., Rho, O. and Hong, S. (2003) Cures for the Shock Instability: Development of a Shock-Stable Roe Scheme. Journal of Computational Physics, 185, 342-374.
[17]	Lax, P.D. (1954) Weak Solutions of Nonlinear Hyperbolic Equations and Their Numerical Computation. Communications on Pure and Applied Mathematics, 7, 159-193. https://doi.org/10.1002/cpa.3160070112
[18]	Toro, E.F. (2009) Riemann Solver and Numerical Methods for Fluid Dynamics. 3rd edition, Springer-Verlag, New York.
[19]	Hu, F. (2017) Conservative and Easily Implemented Finite Volume Semi-Lagran gian WENO Methods for 1D and 2D Hyperbolic Conservation Laws. J. Applied Mathematics and Physics, 5, 59-82. https://doi.org/10.4236/jamp.2017.51008
[20]	Woodward, P. and Colella, P. (1984) The Numerical Simulation of Two-Dimen sional Fluid Flow with Strong Shocks. Journal of Computational Physics, 54, 115-173.