A new type of third-order finite volume multi-resolution WENO schemes on tetrahedral meshes

Journal Pre-proof A new type of third-order finite volume multi-resolution WENO schemes on tetrahedral meshes

Jun Zhu, Chi-Wang Shu

PII:

S0021-9991(19)30917-9

DOI:

https://doi.org/10.1016/j.jcp.2019.109212

Reference:

YJCPH 109212

To appear in:

Journal of Computational Physics

Received date:

9 September 2019

Revised date:

16 December 2019

Accepted date:

17 December 2019

Please cite this article as: J. Zhu, C.-W. Shu, A new type of third-order finite volume multi-resolution WENO schemes on tetrahedral meshes, J. Comput. Phys. (2020), 109212, doi: https://doi.org/10.1016/j.jcp.2019.109212.

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier.

Highlights • A new finite volume multi-resolution WENO scheme on 3D tetrahedral meshes is designed. • Advantages include small number of stencils, arbitrary linear weights provided that their sum is one, and good steady state convergence properties. • Numerical examples are provided to demonstrate the good performance of this new scheme.

A new type of third-order finite volume multi-resolution WENO schemes on tetrahedral meshes Jun Zhu1 and Chi-Wang Shu2 Abstract In this continuing paper of [J. Comput. Phys., 375 (2018), 659-683; J. Comput. Phys., 392 (2019), 19-33], we design a new third-order finite volume multi-resolution weighted essentially non-oscillatory (WENO) scheme for solving hyperbolic conservation laws on tetrahedral meshes. We only use the information defined on a hierarchy of nested central spatial stencils without introducing any equivalent multi-resolution representation. Comparing with classical third-order finite volume WENO schemes [Commun. Comput. Phys., 5 (2009), 836-848] on tetrahedral meshes, the crucial advantages of such new multi-resolution WENO schemes are their simplicity and compactness with the application of only three unequal-sized central stencils for reconstructing unequal degree polynomials in the WENO type spatial procedures, their easy choice of arbitrary positive linear weights without considering the topology of the tetrahedral meshes, their optimal order of accuracy in smooth regions, and their suppression of spurious oscillations near strong discontinuities. The linear weights of such new multi-resolution WENO scheme can be any positive numbers on the condition that their sum is one. This is the first time that a series of unequal-sized hierarchical central spatial stencils are used in designing high-order finite volume WENO scheme on tetrahedral meshes. By performing such new spatial reconstruction procedures and adopting a third-order TVD Runge-Kutta method for time discretization, the occupied memory is decreased and the computing efficiency is increased. This new third-order finite volume multi-resolution WENO scheme is suitable for large scale engineering applications and could maintain good convergence property for steady-state problems on tetrahedral meshes. Benchmark examples are computed to demonstrate the robustness and good performance of these new finite volume WENO schemes. Key Words: multi-resolution WENO scheme, unequal-sized central stencil, finite volume, large scale engineering application, steady-state problems. AMS(MOS) subject classification: 65M08, 35L65 1

College of Science, Nanjing University of Aeronautics and Astronautics, Nanjing, Jiangsu 210016, P.R. China. E-mail: [email protected]. Research was supported by Science Challenge Project, No. TZ2016002 and NSFC grant 11872210. 2 Division of Applied Mathematics, Brown University, Providence, RI 02912, USA. E-mail: chiwang [email protected]. Research was supported by NSF grant DMS-1719410.

1

1

Introduction

In this paper, we design a new third-order finite volume multi-resolution weighted essentially non-oscillatory (WENO) scheme to solve three-dimensional hyperbolic conservation laws  ut + f (u)x + g(u)y + r(u)z = 0, (1.1) u(x, y, z, 0) = u0 (x, y, z), on tetrahedral meshes. Let us first mention some advantages of this new multi-resolution WENO scheme, which is an extension of [62, 63] from one-dimensional and two-dimensional structured or unstructured meshes to three-dimensional tetrahedral meshes. The first is that the linear weights can be arbitrarily set as any positive numbers with one requirement that their sum equals to one. The second is its simplicity and easy extension to unstructured meshes. The third is that a multi-resolution style hierarchy of nested central spatial stencils is used in the WENO spatial reconstruction, which is different from other classical ENO and WENO schemes on non-uniform meshes [1, 2, 11, 16, 20, 24, 33, 34, 41, 44, 45, 50, 51, 52, 55, 56, 61]. The fourth is that the number of spatial stencils is three, which is much smaller than sixteen of the classical third-order finite volume WENO scheme proposed by Zhang and Shu [56] on tetrahedral meshes. The last but not the least feature is that this thirdorder finite volume multi-resolution WENO scheme can obtain optimal third-order accuracy in smooth regions and could degrade gradually to second-order or ultimately to first-order accuracy in non-smooth regions in a non-oscillatory fashion. The performance of this finite volume multi-resolution WENO scheme is good for some benchmark large scale engineering applications and steady-state problems on tetrahedral meshes. In recent decades, essentially non-oscillatory (ENO) and WENO schemes are popular high-order numerical schemes designed for solving hyperbolic conservation laws (1.1) with various wave structures on structured or unstructured meshes in multi-dimensions. We first briefly review the history of these developments. In 1987, Harten et al. designed finite volume ENO schemes [31, 32]. In 1988 and 1989, Shu and Osher developed an efficient implementation of finite difference ENO schemes with total variation diminishing (TVD) time discretizations [47, 48]. In 1994, Liu et al. [40] presented a third-order finite volume WENO scheme. In 1996, Jiang and Shu [35] designed third-order and fifth-order finite difference WENO schemes in multi-dimension with a general framework for designing new smoothness indicators and nonlinear weights. Some classical finite volume WENO schemes were designed [24, 34, 38, 42, 45] on different kinds of computational meshes by adopting above mentioned reconstruction methodologies. These classical ENO and WENO schemes use the idea of adaptive stencils to automatically achieve high-order accuracy in smooth regions and keep essentially non-oscillatory property in non-smooth regions. The crucial merit of the WENO 2

schemes superior to the ENO schemes is that they give a convex combination relationship of all equal degree polynomials defined on different equal-sized spatial stencils and serve as a building block for higher-order spatial reconstructions. Although ENO and WENO schemes [1, 12, 13, 18, 24, 25] are widely used in many engineering fields, it is difficult to extend such WENO schemes from structured meshes to unstructured meshes, such as triangular meshes or tetrahedral meshes. In 1999, Hu and Shu [34] proposed third-order and fourth-order finite volume WENO schemes on triangular meshes. They gave a new way of computing two-dimensional smoothness indicators which was different to [1, 24] and gave a new way of computing the optimal linear weights at different Gaussian quadrature points on the boundaries of the target triangle. But such procedures of obtaining the optimal linear weights for high-order finite volume WENO schemes are complicated and not easy to be accomplished. This difficulty restricts large scale engineering applications on triangular meshes or tetrahedral meshes. In [56], Zhang and Shu classified two types of WENO schemes designed in the literature: the first type includes some WENO schemes whose order of accuracy is not higher than the reconstruction defined on each small spatial stencil. For this type of WENO schemes, the linear weights could be chosen arbitrarily, and the nonlinear weights are only used to avoid spurious oscillations and to keep stability in non-smooth regions. However, these so called robust WENO schemes do not contribute to an increased accuracy from that of each of the small stencils. Since the linear weights could be set as any arbitrary positive constants with one requirement that their sum equals to one, such robust but low-order WENO schemes [20, 21, 24, 51] are easier to construct and suitable for large scale engineering applications. The second type includes some classical WENO schemes [34, 45, 56], whose order of accuracy is higher than that of the reconstruction defined on each of the small spatial stencils. For example, the thirdorder and fourth-order WENO schemes [34] are based on nine linear polynomials defined on three-cell spatial stencils, and on six quadratic polynomials defined on six-cell spatial stencils, respectively. These high-order WENO schemes are difficult to construct, because their optimal linear weights must be computed differently at different Gaussian quadrature points on the boundaries of the target cell, and these linear weights may become negative [45] or even fail to exist. But they have a more compact spatial stencil than the former type of WENO schemes to obtain the same order of accuracy. If the computational meshes distort greatly or change in time (such as for the space-time adaptive meshes [19, 22, 54] or Arbitrary-Lagrange-Eulerian (ALE) problems [8]), we would need to deal with the costly computation of the linear weights and the possible occurrence of the negative linear weights, possibly at every time step. Liu and Zhang [41] combined two types of WENO reconstructions to design robust WENO schemes when facing distorted 3

local mesh geometries or degenerate cases. If the linear weights were significantly negative or non-existent in such circumstances, it would significantly affect the performance of the finite volume WENO schemes. When there are no linear weights at some specific quadrature points on the boundaries of the target cell, the WENO reconstruction procedure will fail to increase the order of accuracy to the optimal one. Such problem of non-existent linear weights for third-order WENO reconstruction at the center of the target cell was addressed in the central/compact WENO (CWENO) schemes [6, 9, 11, 16, 36, 38, 39, 44]. For overcoming the drawback of WENO spatial reconstruction that the linear weights may not exist, Zhu and Qiu have also developed new fifth-order finite difference or finite volume WENO schemes [58, 59, 60, 61] with a series of unequal-sized spatial stencils. In recent decades, the multi-resolution method was proposed by Harten in a series of papers [26, 27, 28, 29, 30], for reducing the computational costs of high-resolution schemes. Since numerical solution of hyperbolic conservation laws might contain strong discontinuities in small and isolated regions, and might be smooth in the remaining large regions, the multi-resolution method could focus its effort in non-smooth regions which contain shocks or contact discontinuities. So Harten modified the original multi-resolution method to be a new version [26, 27, 28, 29, 30] for solving the hyperbolic equations. After that, Abgrall [3, 4] and with Harten [5] designed multi-resolution methods on unstructured meshes. In 2001, Dahmen et al. [17] analyzed the multi-resolution methods for the conservation laws. In 2003, Chiavassa et al. [15] proposed the multi-resolution-based adaptive schemes for simulating the hyperbolic conservation laws. In 2007, B¨ urger et al. [10] proposed a new fifth-order WENO scheme with a multi-resolution method for multi-species kinematic flow models. All in all, the main purpose of applying the multi-resolution method is to focus the computational effort mainly in the small non-smooth regions containing strong shocks or contact discontinuities. Following the original idea of classical CWENO schemes [38, 39] and the new finite difference or finite volume multi-resolution WENO schemes [62, 63], we extend such new multiresolution WENO schemes from one-dimensional and two-dimensional structured or unstructured meshes to three-dimensional tetrahedral meshes, without introducing any equivalent multi-resolution representation as in [26, 27, 28, 29, 30]. The basic flowcharts of the spatial reconstruction are briefly narrated as follows. We select a big central spatial stencil which contains no fewer than ten tetrahedral cells including the target cell, and then reconstruct a quadratic polynomial based on the information of the conservative variables defined on each tetrahedral cell. Hereafter, we select a smaller central spatial stencil which contains no fewer than four tetrahedral cells including the target cell, and reconstruct a linear polynomial based on the information of the conservative variables defined on these tetrahedral cells. 4

Then we select a smallest central spatial stencil which only contains the target cell, and reconstruct a zeroth degree polynomial based on it. After that, the quadratic polynomial needs to be modified, so as to keep third-order approximation in smooth regions [34, 46, 56] and sustain essentially non-oscillatory property in non-smooth regions for the eventual WENO convex combination. After performing these modifications, any positive linear weights could be chosen provided their summation is one. After the computation of the smoothness indicators, and the applications of a new formulation of the nonlinear weights and a third-order TVD Runge-Kutta time discretization method [47], a new third-order finite volume multiresolution WENO scheme is obtained both in space and time. Moreover, the number of the unequal-sized central spatial stencils is two for the second-order scheme and three for the third-order scheme, which is much smaller than that of the same order classical WENO scheme [56] on tetrahedral meshes. The organization of this paper is as follows. The methodology of the new third-order finite volume multi-resolution WENO scheme is described in Section 2. Extensive benchmark one-dimensional, two-dimensional, and three-dimensional numerical examples are provided in Section 3 for verifying the stability, convergence of steady-state problems, and numerical accuracy of the new WENO scheme. Finally, concluding remarks are given in Section 4.

2

Finite volume multi-resolution WENO scheme

We study three-dimensional conservation laws (1.1) on tetrahedral meshes and integrate (1.1) over the target cell 0 to get a semi-discrete finite volume formulation  d¯ u0 (t) 1 + F · nds = 0, (2.1) dt |0 | ∂0  where u¯0 (t) = |10 | 0 u(x, y, z, t)dxdydz, F = (f, g, r), ∂0 is the boundary of the target cell 0 , |0 | is the volume of the target cell 0 , and n = (nx , ny , nz )T denotes the outward unit normal to the boundary of the target cell Δ0 . The surface integrals in (2.1) are discretized by a six-point quadrature integration formula on every triangular face (if we define three vertexes (x1 , y1 , z1 ), (x2 , y2 , z2 ) and (x3 , y3 , z3 ) for a triangular element, the associated six-point quadrature points are (xG1 , yG1 , zG1 ) = (λ1 x1 +λ2 (x2 + x3 ), λ1 y1 +λ2 (y2 + y3 ), λ1 z1 +λ2 (z2 + z3 )), (xG2 , yG2 , zG2 ) = (λ2 (x1 + x3 )+λ1 x2 , λ2 (y1 + y3 )+λ1 y2 , λ2 (z1 + z3 )+λ1 z2 ), (xG3 , yG3 , zG3 ) = (λ2 (x1 + x2 )+λ1 x3 , λ2 (y1 + y2 )+λ1 y3 , λ2 (z1 + z2 )+λ1 z3 ), (xG4 , yG4 , zG4 ) = (β1 x1 +β2 (x2 + x3 ), β1 y1 +β2 (y2 + y3 ), β1 z1 +β2 (z2 + z3 )), (xG5 , yG5 , zG5 ) = (β2 (x1 + x3 )+β1 x2 , β2 (y1 + y3 )+β1 y2 , β2 (z1 + z3 )+β1 z2 ), and (xG6 , yG6 , zG6 ) = (β2 (x1 + x2 )+β1 x3 , β2 (y1 + y2 )+β1 y3 , β2 (z1 + z2 )+β1 z3 ), where λ1 = 0.816847572980459, λ2 = 0.091576213509771, β1 = 0.108103018168070, β2 = 0.445948490915965, and the quadrature weights are σ = 5

0.109951743655322,  = 1, 2, 3 and σ = 0.223381589678011,  = 4, 5, 6, respectively):  ∂0

F · nds ≈

4 

|∂0 |

=1

6  =1

σ F (u(xG , yG , zG , t)) · n .

(2.2)

F (u(xG , yG , zG , t)) · n ,  = 1, ..., 6,  = 1, 2, 3, 4 are reformulated by numerical fluxes such as the simple Lax-Friedrichs flux F (u(xG , yG , zG , t)) · n ≈ 1 [(F (u+ (xG , yG , zG , t)) + F (u− (xG , yG , zG , t))) · n − (2.3) 2 α(u+ (xG , yG , zG , t) − u− (xG , yG , zG , t))],  = 1, ..., 6,  = 1, 2, 3, 4. Here α is an upper bound for the eigenvalues of the Jacobian in the n direction, and u− and u+ are the values of u inside and outside the boundaries of the target tetrahedral cell (inside the neighboring tetrahedral cell) at different quadrature points and |∂0 |,  = 1, 2, 3, 4 are the areas of the triangular elements. We give the flowchart of the new third-order finite volume multi-resolution WENO scheme on tetrahedral meshes and simply omit the variable t in the following, if it does not cause confusion. The reconstruction of the function u(x, y, z) in the target cell 0 , which is used to compute the approximation at different quadrature points (xG , yG , zG ),  = 1, ..., 6,  = 1, 2, 3, 4 on the boundaries of the target cell 0 , is narrated as follows. Step 1. Select a series of central spatial stencils and reconstruct polynomials of different degrees. Step 1.1. For a second-order spatial approximation, we choose two central spatial stencils

T1 = {0 } and T2 = {0 , 1 , 2 , 3 , 4 }, which includes the target cell Δ0 and its four neighboring tetrahedrons Δ1 , Δ2 , Δ3 , and Δ4 , respectively. It is easy to reconstruct a zeroth degree polynomial q1 (x, y, z) which has the same cell average of u on the target cell 0 , and a linear polynomial q2 (x, y, z) ∈ span{1, x−x01 , y−y01 , z−z01 }, where (x0 , y0 , z0 ) is the |0 | 3

|0 | 3

|0 | 3

barycenter of the target cell 0 , which has the same cell average of u on the target cell 0 and matches the cell averages of u on the other tetrahedrons in the set T2 \ {0 } in a least square sense [34]:  1 q1 (x, y, z)dxdydz = u¯0 , (2.4) |0 | 0 1 |0 |

 0

q2 (x, y, z)dxdydz = u¯0 ,

6

(2.5)

and 2  1  p(x, y, z)dxdydz − u¯ , A = {1, 2, 3, 4}, (2.6) q2 (x, y, z) = argminp | |  ∈A where the minimum is taken over all linear polynomials p satisfying the condition (2.5). Step 1.2. For a third-order spatial approximation, we use a central spatial stencil T3 = {0 , 1 , 2 , 3 , 4 , 11 , 12 , 13 , 21 , 22 , 23 , 31 , 32 , 33 , 41 , 42 , 43 }. This big central spatial stencil T3 includes the target cell Δ0 , its four neighboring tetrahedrons Δ1 , Δ2 , Δ3 , and Δ4 , and the neighboring tetrahedrons of these four tetrahedrons (two layers of neighboring tetrahedrons). In T3 , the subscripts are self-evident, for example, Δ11 , Δ12 , and Δ13 are the three neighboring tetrahedrons of Δ1 other than the target cell Δ0 . Therefore, we can construct a quadratic polynomial q3 (x, y, z) ∈ span{1, x−x01 , y−y01 , z−z01 , (x−x0 )2

(x−x0 )(y−y0 )

(x−x0 )(z−z0 )

(y−y0 )2

(y−y0 )(z−z0 )

(z−z0 )2

|0 | 3

|0 | 3

|0 | 3

, , , , 2 2 2 , 2 2 } on T3 to obtain a third-order |0 | 3 |0 | 3 |0 | 3 |0 | 3 |0 | 3 approximation of the conservative variable u. This quadratic polynomial has the same cell average of u on the target cell 0 and matches the cell averages of u on the other tetrahedrons in the set T3 \ {0 } in a least square sense [34]:   1 1 q3 (x, y, z)dxdydz = u(x, y, z)dxdydz = u¯0 , (2.7) |0 | 0 |0 | 0 2

|0 | 3

and 2  1  p(x, y, z)dxdydz − u¯ , q3 (x, y, z) = argminp |  |  ∈A

(2.8)

A = {1, 2, 3, 4, 11, 12, 13, 21, 22, 23, 31, 32, 33, 41, 42, 43}, where the minimum is taken over all quadratic polynomials p satisfying the condition (2.7). For the purpose of obtaining ten degrees of freedom for q3 (x, y, z), we would like to have a central spatial stencil including at least ten distinct tetrahedrons. However, the quality of the tetrahedral meshes might be poor and some tetrahedrons in T3 might coincide with each other. This might result in the lack of enough (ten) distinct tetrahedrons for reconstructing q3 (x, y, z). In such circumstance, we might want to add some or all neighboring tetrahedrons in the next layer around T3 so as to obtain enough tetrahedrons to obtain the correct accuracy. However, we have not encountered such extreme situation in our numerical experiments. Step 2. Obtain equivalent expressions for these reconstruction polynomials of different degrees. To keep consistent notation, we will denote p1 (x, y, z) = q1 (x, y, z) and obtain p2 (x, y, z) =

1 γ2 ,2

q2 (x, y, z) −

 2 −1 =1

7

γ,2 p (x, y, z), γ2 ,2

(2.9)

2 with =1 γ,2 = 1 and γ2 ,2 = 0 for 2 = 2, 3. We explain them in detail in the following. Step 2.1. For the second-order approximation, with similar ideas for designing CWENO schemes [11, 38, 39], a polynomial p2 (x, y, z) is defined through p2 (x, y, z) =

1 γ1,2 q2 (x, y, z) − p1 (x, y, z), γ2,2 γ2,2

(2.10)

with γ1,2 + γ2,2 = 1 and γ2,2 = 0. Step 2.2. For the third-order approximation, a polynomial p3 (x, y, z) is defined through  γ,3 1 q3 (x, y, z) − p (x, y, z), γ3,3 γ3,3 =1 2

p3 (x, y, z) =

(2.11)

 with 3=1 γ,3 = 1 and γ3,3 = 0. In these expressions, γ,2 for  = 1, ..., 2 and 2 = 2, 3 are the linear weights, respectively. Based on a balance between the sharp and essentially non-oscillatory shock transitions in non-smooth regions and numerical accuracy in smooth regions, following the practice in γ ¯ [20, 57, 58, 59, 64], we set the linear weights as γ,2 = 2,2 , in which γ¯,2 = 10−1 for ¯l,2 l=1 γ  = 1, ..., 2 and 2 = 2, 3. For example, we take γ¯1,2 = 1 and γ¯2,2 = 10 for the second-order approximation; γ¯1,3 = 1, γ¯2,3 = 10, and γ¯3,3 = 100 for the third-order approximation. Step 3. Compute the smoothness indicators β2 , which measure how smooth the functions p2 (x, y, z) for 2 = 2, 3 are in the target cell Δ0 . We use the same recipe for the smoothness indicators as in [34, 35, 46]: β =

r  |l|=1

|0 |

2|l| −1 3



 0

∂ |l| p (x, y, z) ∂xl1 ∂y l2 ∂z l3

2 dxdydz,  = 2, 3,

(2.12)

where l = (l1 , l2 , l3 ), |l| = l1 + l2 + l3 . For  = 2, r equals to 1; for  = 3, r equals to 2. The direct application of (2.12) results in β1 = 0. This would still yield optimal order of accuracy in smooth regions but it would affect the maintenance of sharp shock transitions, as the first order building block is given too much weight. The correct measurement of the smoothness of the target cell would be the smallest (most smooth in some weighted sense) among that of the four first order polynomials containing this target cell, that is, we should magnify the quantity of β1 from zero to a suitable bigger value defined in the following. We construct four polynomials p1, (x, y, z) ∈ span{ x−x1 , | | 3

y−y

1

| | 3

,

z−z

1

| | 3

},  = 1, 2, 3, 4,

u11 − u¯1 , p1,1 (x12 , y12 , z12 )=¯ u12 − u¯1 , p1,1 (x13 , y13 , z13 )=¯ u13 − u¯1 ; satisfying p1,1 (x11 , y11 , z11 )=¯ p1,2 (x21 , y21 , z21 )=¯ u21 − u¯2 , p1,2 (x22 , y22 , z22 )=¯ u22 − u¯2 , p1,2 (x23 , y23 , z23 )=¯ u23 − u¯2 ; p1,3 (x31 , y31 , u31 −¯ u3 , p1,3 (x32 , y32 , z32 )=¯ u32 −¯ u3 , p1,3 (x33 , y33 , z33 )=¯ u33 −¯ u3 ; p1,4 (x41 , y41 , z41 )=¯ u41 −¯ u4 , z31 )=¯ p1,4 (x42 , y42 , z42 )=¯ u42 − u¯4 , p1,4 (x43 , y43 , z43 )=¯ u43 − u¯4 . Here (x , y , z ) are the barycenters of 8

Δ ,  = 1, 11, 12, 13; 2, 21, 22, 23; 3, 31, 32, 33; 4, 41, 42, 43, respectively. Then we apply (2.12) ¯ 1,1 =λ ¯ 1,2 =λ ¯ 1,3 =λ ¯ 1,4 =1. to obtain associated smoothness indicators β1, ,  = 1, 2, 3, 4. We set λ λ1,1 = λ1,3 ,

¯ 1,1 λ ¯ 1,2 +λ ¯ 1,3 +λ ¯ 1,4 , ¯ 1,1 +λ λ



σ = λ1, 1 +

(

λ1,2 =

¯ 1,2 λ ¯ 1,2 +λ ¯ 1,3 +λ ¯ 1,4 , ¯ 1,1 +λ λ

λ1,3 =

¯ 1,3 λ ¯ 1,2 +λ ¯ 1,3 +λ ¯ 1,4 , ¯ 1,1 +λ λ

|β1,1 −β1,2 |+|β1,1 −β1,3 |+|β1,1 −β1,4 |+|β1,2 −β1,3 |+|β1,2 −β1,4 |+|β1,3 −β1,4 | 2 ) 6

β1, +ε

λ1,4 = 1-λ1,1 -λ1,2 -

 ,  = 1, 2, 3, 4, and

σ = σ1 + σ2 + σ3 + σ4 , where ε = 10−6 . Then we obtain ⎛   4 2 ⎞ |l|   σ ∂ ⎝|0 ||l| ⎠. p1, (x, y, z) β1 = l l l 1 2 3 ∂x ∂y ∂z σ =1

(2.13)

|l|=1

Step 4. Compute the nonlinear weights based on the linear weights and the smoothness indicators. We adopt the WENO-Z recipe as shown in [7, 14], with τ2 for 2 = 2, 3 defined as related to the absolute difference between the smoothness indicators: 2  2 −1 |β − β |   2 =1 τ2 = , 2 = 2, 3. (2.14) 2 − 1 The nonlinear weights are given as ω ¯ 1 ,2

ω1 ,2 = 2

=1

ω ¯ ,2



, ω ¯ 1 ,2 = γ1 ,2

τ2 1+ ε + β1

 , 1 = 1, ..., 2 ; 2 = 2, 3.

(2.15)

Here ε = 10−6 is chosen in the numerical tests. Step 5. The new final reconstruction polynomial for the value u(x, y, z) is given by Q(x, y, z) =

2 

ω,2 p (x, y, z), 2 = 2, 3,

(2.16)

=1

for the second-order and third-order approximations, respectively. The approximations can then be given as u− (xG , yG , zG ) ≈ Q(xG , yG , zG ),  = 1, ..., 6,  = 1, 2, 3, 4,

(2.17)

at different quadrature points on the boundaries of the target cell 0 for the second-order and third-order approximations, respectively. Step 6. A third-order TVD Runge-Kutta time discretization method [47] ⎧ (1) ⎨ u = un + tL(un ), (2.18) u(2) = 3 un + 1 u(1) + 1 tL(u(1) ), ⎩ n+1 4 1 n 4 2 (2) 4 2 (2) u = 3 u + 3 u + 3 tL(u ), is used to solve (2.1). Finally, the fully discrete scheme both in space and time is designed on tetrahedral meshes. 9

3

Numerical results

In this section, we simulate some benchmark numerical examples to testify the performance of the new third-order finite volume multi-resolution WENO scheme with only three unequal-sized central spatial stencils on tetrahedral meshes. For systems of the compressible Euler equations, the spatial reconstructions are used in the local characteristic directions to avoid spurious oscillations near strong discontinuities. The CFL number is set as 0.6. For the temporal discretization, the third-order TVD Runge-Kutta time discretization method [47] (2.18) is adopted for all examples. First, for evaluating whether arbitrary choice of the positive linear weights would always sustain optimal third-order accuracy in smooth regions, γ ¯

we set four different types of the linear weights as γ,3 = 3 ,3γ¯ , in which (1) γ¯1,3 = 1, l=1 l,3 γ¯2,3 = 10, and γ¯3,3 = 100; (2) γ¯1,3 = 1, γ¯2,3 = 10, and γ¯3,3 = 10; (3) γ¯1,3 = 1, γ¯2,3 = 1, and γ¯3,3 = 100; (4) γ¯1,3 = 1, γ¯2,3 = 1, and γ¯3,3 = 1 in all numerical accuracy tests. Then, we will use only γ¯1,3 = 1, γ¯2,3 = 10, and γ¯3,3 = 100 in the other test cases, unless specified otherwise. Example 3.1. We solve the following three-dimensional linear scalar equation ut + ux + uy + uz = 0,

(3.1)

with the computational domain [−2, 2] × [−2, 2] × [−2, 2] on uniform tetrahedral meshes. The initial condition is u(x, y, z, 0) = sin(π(x + y + z)/2) and periodic boundary conditions are applied in different directions. The final time is t = 1. The errors and numerical orders of accuracy for the third-order finite volume multi-resolution WENO scheme with different types of the linear weights are shown in Table 3.1. We can see that the new third-order WENO scheme keeps the designed order of accuracy. We do notice that the magnitude of the absolute L1 and L∞ errors are modestly different for different types of the linear weights. Example 3.2. We solve the following three-dimensional Burgers’ equation  2  2  2 u u u + + = 0, ut + 2 x 2 y 2 z

(3.2)

with the computational domain [−3, 3] × [−3, 3] × [−3, 3] on the uniform tetrahedral meshes. The initial condition is u(x, y, z, 0) = 0.5 + sin(π(x + y + z)/3) and periodic boundary conditions are applied in different directions. The final time is t = 0.5/π 2 . The errors and numerical orders of accuracy for the new multi-resolution WENO scheme with four different types of the linear weights are shown in Table 3.2. We can also see that the new third-order scheme keeps the designed order of accuracy, and the magnitude of the absolute truncation errors are modestly different for different types of the linear weights on the same mesh levels once again. 10

Table 3.1: ut + ux + uy + uz = 0. u(x, y, z, 0) = sin(π(x + y + z)/2). Periodic boundary conditions in each direction. T = 1. L1 and L∞ errors. Uniform tetrahedral mesh. tetrahedrons 750 6000 20250 48000 93750 tetrahedrons 750 6000 20250 48000 93750

WENO scheme (1) L1 error order L∞ error 3.09E-1 5.04E-1 4.77E-2 2.70 9.35E-2 1.36E-2 3.08 2.80E-2 5.83E-3 2.97 9.88E-3 3.01E-3 2.96 4.89E-3 WENO scheme (3) L1 error order L∞ error 3.16E-1 5.11E-1 5.12E-2 2.63 9.62E-2 1.50E-2 3.01 2.86E-2 6.51E-3 2.92 9.95E-3 3.39E-3 2.92 4.88E-3

order 2.43 2.97 3.62 3.15 order 2.41 2.99 3.67 3.19

WENO scheme (2) L1 error order L∞ error 3.58E-1 5.82E-1 5.98E-2 2.58 1.35E-1 1.38E-2 3.61 4.29E-2 5.74E-3 3.06 1.60E-2 3.00E-3 2.91 6.64E-3 WENO scheme (4) L1 error order L∞ error 5.12E-1 8.22E-1 2.65E-1 0.95 4.57E-1 8.39E-2 2.84 1.88E-1 2.98E-2 3.60 7.65E-2 1.55E-2 2.92 3.32E-2

Example 3.3. We solve the three-dimensional Euler equations ⎞ ⎞ ⎞ ⎛ ⎛ ⎛ ⎛ ρ ρu ρv ρw ⎟ ⎜ ρu ⎟ ⎜ ρu2 + p ⎟ ⎜ ⎜ ρuv ρuw ⎟ ⎟ ⎟ ⎜ ⎜ ⎜ ∂ ⎜ ⎟ + ∂ ⎜ ρv 2 + p ⎟ + ∂ ⎜ ⎜ ρv ⎟ + ∂ ⎜ ρvw ρvu ⎟ ⎟ ⎟ ⎜ ⎜ ⎜ ∂t ⎜ ⎝ ρw ⎠ ∂x ⎝ ρwu ⎠ ∂y ⎝ ρwv ⎠ ∂z ⎝ ρw2 + p E u(E + p) v(E + p) w(E + p)

order 2.10 2.84 3.43 3.94 order 0.85 2.19 3.13 3.74

⎞ ⎟ ⎟ ⎟ = 0, ⎟ ⎠

(3.3)

in which ρ is density, u is x-direction velocity, v is y-direction velocity, w is z-direction velocity, respectively, E is total energy, and p is pressure, related to the conserved variables p by the equation of state for ideal gas E = γ−1 + 12 (u2 + v 2 + w2 ) with γ = 1.4. The initial conditions are ρ(x, y, z, 0) = 1 + 0.2 sin(π(x + y + z)/3), u(x, y, z, 0) = 1, v(x, y, z, 0) = 1, w(x, y, z, 0) = 1, and p(x, y, 0) = 1 with the computational domain [−3, 3] × [−3, 3] × [−3, 3] on the uniform tetrahedral meshes. Periodic boundary conditions are applied in all three directions. The final time is t = 1. The errors and numerical orders of accuracy of the density for the third-order finite volume multi-resolution WENO scheme with different types of the linear weights are shown in Table 3.3. Similar to the previous scalar examples, we can see that the new scheme keeps its optimal third-order accuracy in this benchmark example again. Example 3.4. We solve the three-dimensional Euler equations (3.3) with the Riemann initial condition for the Lax problem [37]: (ρ, u, v, w, p, γ)T =



(0.445, 0.698, 0, 0, 3.528, 1.4)T , (0.5, 0, 0, 0, 0.571, 1.4)T ,

(x, y, z)T ∈ [−0.5, 0) × [−0.02, 0.02] × [−0.02, 0.02], (x, y, z)T ∈ [0, 0.5] × [−0.02, 0.02] × [−0.02, 0.02].

11

(3.4)

Table 3.2: ut +

 2 u 2

x

+

 2 u 2

y

+

 2 u 2

z

= 0. u(x, y, z, 0) = 0.5 + sin(π(x + y + z)/3). Periodic

boundary conditions in each direction. T = 0.5/π 2 . L1 and L∞ errors. Uniform tetrahedral mesh. tetrahedrons 750 6000 20250 48000 93750 tetrahedrons 750 6000 20250 48000 93750

WENO scheme (1) L1 error order L∞ error 2.42E-2 6.77E-2 3.04E-3 2.99 1.67E-2 7.61E-4 3.42 2.03E-3 3.26E-4 2.94 7.90E-4 1.68E-4 2.97 4.17E-4 WENO scheme (3) L1 error order L∞ error 2.61E-2 7.22E-2 3.22E-3 3.02 1.81E-2 8.39E-4 3.32 2.23E-3 3.63E-4 2.91 8.84E-4 1.87E-4 2.95 4.74E-4

order 2.02 5.20 3.29 2.86 order 1.99 5.17 3.22 2.79

WENO scheme (2) L1 error order L∞ error 2.97E-2 7.14E-2 4.55E-3 2.71 2.26E-2 1.18E-3 3.32 7.15E-3 4.63E-4 3.26 2.24E-3 1.68E-4 4.53 1.04E-3 WENO scheme (4) L1 error order L∞ error 6.31E-2 1.19E-1 2.01E-2 1.65 5.48E-2 6.33E-3 2.85 2.48E-2 2.18E-3 3.71 9.88E-3 9.72E-4 3.62 4.11E-3

order 1.66 2.84 4.03 3.43 order 1.13 1.95 3.20 3.93

The three-dimensional tetrahedral finite volume multi-resolution WENO scheme is applied to the one-dimensional shock tube problem. The solution lies in the domain of [−0.5, 0.5] × [−0.02, 0.02] × [−0.02, 0.02] with a tetrahedralization of 101 vertices in the x-direction and 5 vertices in the y-direction and z-direction, respectively. The velocities in the y-direction and z-direction are set as 0 and periodic boundary conditions are applied in these directions. The final time is t = 0.16. We present the exact solution and the computed density ρ obtained with the new third-order finite volume multi-resolution WENO scheme in Figure 3.1. We observe that the computational results obtained by the new WENO scheme are good. Example 3.5. We solve the three-dimensional Euler equations (3.3) with the Riemann initial condition for the Sod problem [49]: (ρ, u, v, w, p, γ)T =



(1, 0, 0, 0, 2.5, 1.4)T , (0.125, 0, 0, 0, 0.25, 1.4)T ,

(x, y, z)T ∈ [−5, 0) × [−0.2, 0.2] × [−0.2, 0.2], (x, y, z)T ∈ [0, 5] × [−0.2, 0.2] × [−0.2, 0.2].

(3.5)

The solution lies in the domain of [−5, 5]×[−0.2, 0.2]×[−0.2, 0.2] with a tetrahedralization of 101 vertices in the x-direction and 5 vertices in the y-direction and z-direction, respectively. The velocities in the y-direction and z-direction are set as 0 and periodic boundary conditions are applied in these directions. The final time is t = 2. We present the exact solution and the computed density ρ obtained with the new third-order finite volume multi-resolution WENO scheme in Figure 3.2. The numerical results computed by the WENO scheme are good for this one-dimensional test case. 12

1.4 DEN

1.2

Y X

1

Density

1

Density

Z

0.88373 0.462906 0.0420824

1.2

0.8

0.6

0.8

0.6

0.4

0.4

-0.4

-0.2

0

0.2

0.4

0.2

0.4

X

X

Figure 3.1: The Lax problem. T=0.16. From left to right: density cutting-plot along y = z = 0; density zoomed in; density surface cutting-plot along z = 0. Solid line: the exact solution; squares: the results of WENO scheme. The mesh points on the boundary are uniformly distributed with cell length Δx = Δy = Δz = 1/100.

1.2

DEN

1

0.681181 0.356809 0.0324372

1

Z

Y X

0.8

Density

Density

0.8

0.6

0.4

0.4

0.2

0

0.6

0.2

-4

-2

0

2

4

-2

0

X

2

4

X

Figure 3.2: The Sod problem. T=2. From left to right: density cutting-plot along y = z = 0; density zoomed in; density surface cutting-plot along z = 0. Solid line: the exact solution; squares: the results of WENO scheme. The mesh points on the boundary are uniformly distributed with cell length Δx = Δy = Δz = 10/100.

13

Table 3.3: 3D-Euler equations: initial data ρ(x, y, z, 0) = 1 + 0.2 sin(π(x + y + z)/3), u(x, y, z, 0) = 1, v(x, y, z, 0) = 1, w(x, y, z, 0) = 1, and p(x, y, z, 0) = 1. Periodic boundary conditions in each direction. T = 1. L1 and L∞ errors. Uniform tetrahedral mesh. tetrahedrons 750 6000 20250 48000 93750 tetrahedrons 750 6000 20250 48000 93750

WENO scheme (1) L1 error order L∞ error 6.55E-2 1.04E-1 1.31E-2 2.32 2.12E-2 4.17E-3 2.82 6.61E-3 1.79E-3 2.94 2.84E-3 9.28E-4 2.95 1.50E-3 WENO scheme (3) L1 error order L∞ error 6.60E-2 1.05E-1 1.31E-2 2.33 2.11E-2 4.19E-3 2.82 6.62E-3 1.80E-3 2.94 2.84E-3 9.33E-4 2.95 1.50E-3

order 2.30 2.88 2.94 2.86 order 2.32 2.86 2.94 2.86

WENO scheme (2) L1 error order L∞ error 6.99E-2 1.12E-1 1.33E-2 2.39 2.30E-2 4.18E-3 2.87 6.90E-3 1.79E-3 2.94 2.91E-3 9.28E-4 2.95 1.51E-3 WENO scheme (4) L1 error order L∞ error 8.76E-2 1.40E-1 1.95E-2 2.16 3.72E-2 5.07E-3 3.33 9.85E-3 2.11E-3 3.04 3.80E-3 1.09E-3 2.95 1.82E-3

order 2.28 2.97 2.99 2.95 order 1.91 3.28 3.31 3.30

Example 3.6. A higher order scheme would show its advantage when the solution contains both shocks and complex smooth region structures. A typical example for this is the problem of shock interaction with entropy waves [46]. We solve the three-dimensional Euler equations (3.3) with a moving Mach number 3 shock interacting with sine waves in density: (ρ, u, v, w, p, γ)T = (3.857143, 2.629369, 0, 0, 10.333333, 1.4)T for (x, y, z)T ∈ [−5, −4) × [−0.05, 0.05]×[−0.05, 0.05]; (ρ, u, v, w, p, γ)T = (1+0.2 sin(5x), 0, 0, 0, 1, 1.4)T for (x, y, z)T ∈ [−4, 5]×[−0.05, 0.05]×[−0.05, 0.05]. The solution lies in the computational domain [−5, 5]× [−0.05, 0.05] × [−0.05, 0.05] with a tetrahedralization of 401 vertices in the x-direction and 5 vertices in the y-direction and z-direction, respectively. The velocities in the y-direction and z-direction are set as 0 and periodic boundary conditions are applied in these directions. The computed density ρ is plotted at t = 1.8 against the referenced “exact” solution which is a converged solution computed by the one-dimensional fifth-order finite difference WENO scheme [35] with 2000 grid points in Figure 3.3. The new third-order finite volume multi-resolution WENO scheme could get good resolution for this benchmark example. Example 3.7. We now consider the interaction of two blast waves [53]. The initial conditions are (ρ, u, v, w, p, γ)

T

⎧ ⎨ (1, 0, 0, 0, 103 , 1.4)T , = (1, 0, 0, 0, 10−2 , 1.4)T , ⎩ (1, 0, 0, 0, 102 , 1.4)T ,

(x, y, z)T ∈ [0, 0.1] × [−0.005, 0.005] × [−0.005, 0.005], (x, y, z)T ∈ (0.1, 0.9] × [−0.005, 0.005] × [−0.005, 0.005], (x, y, z)T ∈ (0.9, 1] × [−0.005, 0.005] × [−0.005, 0.005].

(3.6)

The solution of the Euler equations (3.3) lies in the domain of [0, 1] × [−0.005, 0.005] × 14

DEN

Z Y

3.07343 1.60989 0.146354

4

X

Density

Density

4 3

2

3 1 -4

-2

0

2

4

-2

0

X

2

X

Figure 3.3: The shock density wave interaction problem. T=1.8. From left to right: density cutting-plot along y = z = 0; density zoomed in; density surface cutting-plot along z = 0. Solid line: the exact solution; squares: the results of WENO scheme. The mesh points on the boundary are uniformly distributed with cell length Δx = Δy = Δz = 10/400. DEN

6

6

3.73675 1.95734 0.17794

Z Y

X

5

4

4

Density

Density

5

3

3

2

2

1

1

0

0

0.2

0.4

0.6

0.8

1

0.6

0.8

X

X

Figure 3.4: The blast wave problem. T=0.038. From left to right: density cutting-plot along y = z = 0; density zoomed in; density surface cutting-plot along z = 0. Solid line: the exact solution; squares: the results of WENO scheme. The mesh points on the boundary are uniformly distributed with cell length Δx = Δy = Δz = 1/400. [−0.005, 0.005] with a tetrahedralization of 401 vertices in the x-direction and 5 vertices in the y-direction and z-direction, respectively. The velocities in the y-direction and z-direction are set as 0, and periodic boundary conditions are applied in these directions. The computed density ρ is plotted at t = 0.038 against the reference “exact” solution which is a converged solution computed by the one-dimensional fifth-order finite difference WENO scheme [35] with 2000 grid points in Figure 3.4. The new third-order multi-resolution WENO scheme could get good performance again. Example 3.8. We solve the same three-dimensional nonlinear Burgers’ equation (3.2) with the same initial condition u(x, y, z, 0) = 0.5 + sin(π(x + y + z)/3), except that we plot the results at t = 5/π 2 when a shock has already appeared in the solution. We show the contours 15

Z Y

1.5

X U

2

1

-1

-2 2

Z

0

1.4 1.3 1.2 1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1 -0.2 -0.3 -0.4

1

U

3

0.5

0

1

Y

0 -1 -2 -3-3

-2

-1

0

1

2

3

-3

-0.5 -6

-4

-2

0

2

4

6

X+Y

X

Figure 3.5: Burgers’ equation. T = 5/π 2 . Left: contour plot on the surface; right: 1D cutting-plot along x = y and z = 0 with squares representing the results of WENO scheme and the line is the exact solution. on the surface and one-dimensional cutting-plot along x = y and z = 0 of the solutions by the third-order finite volume multi-resolution WENO scheme in Figure 3.5. We can see that the scheme gives non-oscillatory shock transitions for this scalar problem. Example 3.9. The transonic flow over the Onera M6 wing problem [23]. This problem is a classic CFD validation case for external flows because of its simple geometry combined with complexities of transonic flow. This benchmark test case uses the following flow conditions:  M∞ = 0.84 and angle of attack α = 3.06◦ . The computational domain is x2 + y 2 + z 2 ≤ 16 and z ≥ 0, which consists of 143645 tetrahedrons and 1311 triangles over the surface. The surface mesh and Mach number contours are shown in Figure 3.6. The reduction of density residual as a function of the number of iterations is also shown in the same Figure 3.6. It can also be seen in Figure 3.6 that the average residue of such WENO scheme can settle down to a value around 10−17 , close to machine zero. We can see that the new third-order finite volume multi-resolution scheme performs well in this example. Example 3.10. We use INRIA’s 3D tetrahedral elements for the BTC0 (streamlined body, laminar) test case in project ADIGMA [43] with the initial conditions: M∞ = 0.5 and  angle of attack α = 0◦ . The computational domain is x2 + y 2 + z 2 ≤ 10, which consists of 191753 tetrahedrons and 8244 triangles over the surface. The surface mesh is shown in Figure 3.7. We give Mach number and pressure in Figure 3.7. The reduction of density residual as a function of the number of iterations is also shown in Figure 3.7. It can also be seen in Figure 3.7 that the average residue of such WENO scheme can settle down to a value around 10−16.7 , close to machine zero. It shows that the third-order WENO scheme gives good resolution. Example 3.11. We consider the inviscid Euler transonic flow past a Mig plane whose 16

mach

0.899581 0.864361 X 0.829141 0.793922 Z 0.758702 0.723483 0.688263 0.653044

X Z

-4

Log10 (Residual of Density)

Y

Y

-6

-8

-10

-12

-14

-16

-18

0

20000

40000

60000

Iteration

Figure 3.6: The transonic flow over the Onera M6 wing problem. M∞ = 0.84 and angle of attack α = 3.06◦ . From left to right: Onera M6 wing surface mesh, zoomed in; Mach number contours plot on the surface; the reduction of density residual as a function of the number of iterations. Y

Y

mach

X X Z

Z

1.30686 1.14983 0.992789 0.835752 0.678715 0.521678 0.364641 0.207604

Y

-Cp

X Z

-0.499087 -1.16447 -1.82986

Log10 (Residual of Density)

-4

-6

-8

-10

-12

-14

-16

-18

0

10000

20000

30000

40000

50000

60000

Iteration

Figure 3.7: The transonic flow over the BTC0 problem. M∞ = 0.5, angle of attack α = 0◦ . From left to right and top to bottom: BTC0 surface mesh, zoomed in; Mach number contours plot on the surface; pressure contours plot on the surface; the reduction of density residual as a function of the number of iterations. 17

mach

Z

X

Z

0.811281 0.747361 0.683442 0.619523 0.555603 0.491684

Y

X

Y

-4

0.0734258 -0.293848 -0.661123 -1.0284 -1.39567 -1.76295

Z

X

Y

Log10 (Residual of Density)

-6 -Cp

-8

-10

-12

-14

-16

-18 0

20000

40000

6000

Iteration

Figure 3.8: The transonic flow over the Mig jet plane problem. M∞ = 0.85, angle of attack α = 1◦ . From left to right and top to bottom: Mig jet plane surface mesh, zoomed in; Mach number contours plot on the surface; pressure contours plot on the surface; the reduction of density residual as a function of the number of iterations. surface mesh is also modeled by INRIA. The initial conditions are M∞ = 0.85 and angle of  attack α = 1◦ . The computational domain is x2 + y 2 + z 2 ≤ 50, which consists of 7347 tetrahedrons with 1206 triangles over the surface. The surface mesh is shown in Figure 3.8. Mach number and pressure are shown in Figure 3.8. The reduction of density residual as a function of the number of iterations is also shown in Figure 3.8. It can be seen in Figure 3.8 that the average residue of such WENO scheme can settle down to a value around 10−16.5 , close to machine zero. We can see that the new third-order finite volume multi-resolution WENO scheme performs well for this subsonic example.

4

Concluding remarks

In this paper, a new third-order finite volume multi-resolution WENO scheme is designed for solving the hyperbolic conservation laws on tetrahedral meshes. The main advantages of 18

this new multi-resolution WENO scheme are its simplicity on tetrahedral meshes, compact property, and the application of hierarchical structure with only three unequal-sized spatial stencils in comparison with sixteen equal-sized biased spatial stencils for classical thirdorder finite volume WENO scheme [56]. It is a new way of bypassing the calculation of the optimal linear weights for high-order accuracy and it avoids the performance of dealing with the negative linear weights [45]. This new spatial reconstruction framework of the third-order approximation at any quadrature points by one reconstruction polynomial defined on a big central spatial stencil including at least ten tetrahedrons, one linear polynomial defined on a smaller central stencil including at least four tetrahedrons, and one zeroth degree polynomial defined on the target cell is easily implemented, and it can be generalized straightforwardly to design arbitrary high-order finite volume multi-resolution WENO schemes. Such new multiresolution WENO scheme only depends on the usage of suitable high degree polynomial for obtaining high-order approximation in smooth regions and can degrade gradually to lower order polynomials, ultimately to zeroth degree polynomial, for keeping essentially nonoscillatory property near strong discontinuities. This new third-order finite volume WENO scheme is suitable for solving some benchmark steady-state problems on tetrahedral meshes. The methodology could also be used to design new simple and compact multi-resolution WENO limiters for discontinuous Galerkin methods on unstructured 3D meshes, the study of which is our ongoing work.

References [1] R. Abgrall, On essentially non-oscillatory schemes on unstructured meshes: analysis and implementation, J. Comput. Phys., 114 (1993), 45-58. [2] R. Abgrall, An essentially non-oscillatory reconstruction procedure on finite-element type meshes: application to compressible flows, Comput. Methods Appl. Mech. Eng., 116 (1994), 95-101. [3] R. Abgrall, Multiresolution analysis on unstructured meshes: Application to CFD. In Morton, K. W. (ed.) et al., Numerical methods for fluid dynamics V. Proceedings of the conference, Oxford, UK, April 1995. Oxford: Clarendon Press, 271-277, 1995. [4] R. Abgrall, Multiresolution analysis on unstructured meshes: Application to CFD. Chetverushkin, B. N. (ed.) et al., Experimentation, Modelling and computation in flow, turbulence and combustion. Vol. 2. Proceedings of the 3rd French-Russian-Uzbek workshop on Fluid dynamics, Tashkent, Uzbekistan, April 23-28, 1995. Chichester: Wiley. CMAS: Computational Methods in Applied Sciences, 147-156, 1997. 19

[5] R. Abgrall and A. Harten, Multiresolution representation in unstructured meshes, SIAM J. Numer. Anal., 35(6) (1998), 2128-2146. [6] D.S. Balsara, S. Garain and C.-W. Shu, An efficient class of WENO schemes with adaptive order, J. Comput. Phys., 326 (2016), 780-804. [7] R. Borges, M. Carmona, B. Costa and W.S. Don, An improved weighted essentially non-oscillatory scheme for hyperbolic conservation laws, J. Comput. Phys., 227 (2008), 3191-3211. [8] W. Boscheri and M. Dumbser, A direct Arbitrary-Lagrangian-Eulerian ADER-WENO finite volume scheme on unstructured tetrahedral meshes for conservative and nonconservative hyperbolic systems in 3D, J. Comput. Phys., 275 (2014), 484-523. [9] W. Boscheri, M. Semplice and M. Dumbser, Central WENO subcell finite volume limiters for ADER discontinuous Galerkin schemes on fixed and moving unstructured meshes, Commun. Comput. Phys., 25 (2019), 311-346. [10] R. B¨ urger and A. Kozakevicius, Adaptive multiresolution WENO schemes for multispecies kinematic flow models, J. Comput. Phys., 224 (2007), 1190-1222. [11] G. Capdeville, A central WENO scheme for solving hyperbolic conservation laws on non-uniform meshes, J. Comput. Phys., 227 (2008), 2977-3014. [12] J. Casper, Finite-volume implementation of high-order essentially nonoscillatory schemes in two dimensions, AIAA Journal, 30 (1992), 2829-2835. [13] J. Casper and H.-L. Atkins, A finite-volume high-order ENO scheme for two-dimensional hyperbolic systems, J. Comput. Phys., 106 (1993), 62-76. [14] M. Castro, B. Costa and W.S. Don, High order weighted essentially non-oscillatory WENO-Z schemes for hyperbolic conservation laws, J. Comput. Phys., 230 (2011), 1766-1792. [15] G. Chiavassa, R. Donat and S. M¨ uller, Multiresolution-based adaptive schemes for hyperbolic conservation laws, T. Plewa, T. Linde, V.G. Weiss (Eds.), Adaptive Mesh Refinement-Theory and Applications, Lecture Notes in Computational Science and Engineering, vol. 41, Springer-Verlag, Berlin (2003), 137-159. [16] I. Cravero and M. Semplice, On the accuracy of WENO and CWENO reconstructions of third order on nonuniform meshes, J. Sci. Comput., 67 (2016), 1219-1246. 20

[17] W. Dahmen, B. Gottschlich-M¨ uller and S. M¨ uller, Multiresolution schemes for conservation laws, Numer. Math., 88 (2001), 399-443. [18] Y.L. Du, L. Yuan and Y.H. Wang, A high-order modified finite volume WENO method on 3D Cartesian grids, Commun. Comput. Phys., 26 (2019), 768-784. [19] M. Dumbser, A. Hidalgo and O. Zanotti, High order space-time adaptive ADER-WENO finite volume schemes for non-conservative hyperbolic systems, Comput. Methods Appl. Mech. Eng., 268 (2014), 359-387. [20] M. Dumbser and M. K¨aser, Arbitrary high order non-oscillatory finite volume schemes on unstructured meshes for linear hyperbolic systems, J. Comput. Phys., 221 (2007), 693-723. [21] M. Dumbser, M. K¨aser, V.A. Titarev and E.F. Toro, Quadrature-free non-oscillatory finite volume schemes on unstructured meshes for nonlinear hyperbolic systems, J. Comput. Phys., 226 (2007), 204-243. [22] M. Dumbser, O. Zanotti, A. Hidalgo and D.S. Balsara, ADER-WENO finite volume schemes with space-time adaptive mesh refinement, J. Comput. Phys., 248 (2013), 257286. [23] B. Eisfeld, ONERA M6 wing, in FLOMANIA-A European Initiative on Flow Physics Modelling. Notes on Numerical Fluid Mechanics and Multidisciplinary Design (NNFM), edited by W. Haase, B. Aupoix, U. Bunge, D. Schwamborn. Springer, Berlin, Heidelberg. 94 (2006). [24] O. Friedrichs, Weighted essentially non-oscillatory schemes for the interpolation of mean values on unstructured grids, J. Comput. Phys., 144 (1998), 194-212. [25] A. Harten, Preliminary results on the extension of ENO schemes to two-dimensional problems, in Proceedings, International Conference on Nonlinear Hyperbolic Problems, Saint-Etienne, 1986, Lecture Notes in Mathematics, edited by C. Carasso et al. (Springer-Verlag, Berlin, 1987), 23-40. [26] A. Harten, Multi-resolution analysis for ENO schemes, Institute for Computer Applications in Science and Engineering, NASA Langley Research Center, Hampton, Virginia 23665-5225, Contract No. NAS1-18605, September 1991. [27] A. Harten, Discrete multi-resolution analysis and generalized wavelets, Applied Numerical Mathematics, 12 (1993), 153-192. 21

[28] A. Harten, Adaptive multiresolution schemes for shock computations, Comput. Phys., 115 (1994), 319-338. [29] A. Harten, Multiresolution algorithms for the numerical solution of hyperbolic conservation laws, Communications on Pure and Applied Mathematics, 48 (1995), 1305-1342. [30] A. Harten, Multiresolution representation of data: A general framework, SIAM J. Numer. Anal., 33 (1996), 1205-1256. [31] A. Harten and S.R. Chakravarthy, Multi-dimensional ENO schemes for general geometries, ICASE Report, (1991), 91-76. [32] A. Harten, B. Engquist, S. Osher and S. Chakravarthy, Uniformly high order accurate essentially non-oscillatory schemes, III, J. Comput. Phys., 71 (1987), 231-303. [33] G. Hu, R. Li and T. Tang, A robust WENO type finite volume solver for steady Euler equations on unstructured grids, Commun. Comput. Phys., 9 (2011), 627-648. [34] C. Hu and C.-W. Shu, Weighted essentially non-oscillatory schemes on triangular meshes, J. Comput. Phys., 150 (1999), 97-127. [35] G. Jiang and C.-W. Shu, Efficient implementation of weighted ENO schemes, J. Comput. Phys., 126 (1996), 202-228. [36] O. Kolb, On the full and global accuracy of a compact third order WENO scheme, SIAM J. Numer. Anal., 52 (5) (2014), 2335-2355. [37] P.D. Lax, Weak solutions of nonlinear hyperbolic equations and their numerical computation, Commun. Pure Appl. Math., 7 (1954), 159-193. [38] D. Levy, G. Puppo and G. Russo, Central WENO schemes for hyperbolic systems of conservation laws, Math. Model. Numer. Anal., 33 (1999), 547-571. [39] D. Levy, G. Puppo and G. Russo, Compact central WENO schemes for multidimensional conservation laws, SIAM J. Sci. Comput., 22 (2) (2000), 656-672. [40] X.-D. Liu, S. Osher and T. Chan, Weighted essentially non-oscillatory schemes, J. Comput. Phys., 115 (1994), 200-212. [41] Y. Liu and Y.T. Zhang, A robust reconstruction for unstructured WENO schemes, J. Sci. Comput., 54 (2013), 603-621.

22

[42] J. Qiu and C.-W. Shu, On the construction, comparison, and local characteristic decomposition for high order central WENO schemes, J. Comput. Phys., 183 (2002), 187-209. [43] J. Qiu and J. Zhu, RKDG with WENO Type Limiters, in ADIGMA-A European Initiative on the Development of Adaptive Higher-Order Variational Methods for Aerospace Applications, edited by N. Kroll, H. Bieler, H. Deconinck, V. Couaillier, H. van der Ven, K. S∅ensen, Notes on Numerical Fluid Mechanics and Multidisciplinary Design. Springer, Berlin, Heidelberg. 113 (2010). [44] M. Semplice, A. Coco and G. Russo, Adaptive mesh refinement for hyperbolic systems based on third-order compact WENO reconstruction, J. Sci. Comput., 66 (2016), 692724. [45] J. Shi, C. Hu and C.-W. Shu, A technique of treating negative weights in WENO schemes, J. Comput. Phys., 175 (2002), 108-127. [46] C.-W. Shu, Essentially non-oscillatory and weighted essentially non-oscillatory schemes for hyperbolic conservation laws, in Advanced Numerical Approximation of Nonlinear Hyperbolic Equations, edited by A. Quarteroni, Editor, Lecture Notes in Mathematics, CIME subseries (Springer-Verlag, Berlin/New York); ICASE Report 97-65. [47] C.-W. Shu and S. Osher, Efficient implementation of essentially non-oscillatory shock capturing schemes, J. Comput. Phys., 77 (1988), 439-471. [48] C.-W. Shu and S. Osher, Efficient implementation of essentially non-oscillatory shock capturing schemes II, J. Comput. Phys., 83 (1989), 32-78. [49] G. Sod, A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws, J. Comput. Phys., 27 (1978), 1-31. [50] T. Sonar, On the construction of essentially non-oscillatory finite volume approximations to hyperbolic conservation laws on general triangulations: Polynomial recovery, accuracy and stencil selection, Comput. Methods Appl. Mech. Engrg., 140 (1997), 157181. [51] V.A. Titarev, P. Tsoutsanis, D. Drikakis, WENO schemes for mixed-element unstructured meshes, Commun. Comput. Phys., 8 (2010), 585-609. [52] Z.J. Wang and R.F. Chen, Optimized weighted essentially non-oscillatory schemes for linear waves with discontinuity, J. Comput. Phys., 174 (2001), 381-404.

23

[53] P. Woodward and P. Colella, The numerical simulation of two-dimensional fluid flow with strong shocks, J. Comput. Phys., 54 (1984), 115-173. [54] O. Zanotti and M. Dumbser, A high order special relativistic hydrodynamic and magnetohydrodynamic code with space-time adaptive mesh refinement, Comput. Phys. Commun., 188 (2015), 110-127. [55] Y.T. Zhang and C.-W. Shu, High order WENO schemes for Hamilton-Jacobi equations on triangular meshes, SIAM J. Sci. Comput., 24 (2003), 1005-1030. [56] Y.T. Zhang and C.-W. Shu, Third order WENO scheme on three dimensional tetrahedral meshes, Commun. Comput. Phys., 5 (2009), 836-848. [57] X. Zhong and C.-W. Shu, A simple weighted essentially nonoscillatory limiter for RungeKutta discontinuous Galerkin methods, J. Comput. Phys., 232 (2013), 397-415. [58] J. Zhu and J. Qiu, A new fifth order finite difference WENO scheme for solving hyperbolic conservation laws, J. Comput. Phys., 318 (2016), 110-121. [59] J. Zhu and J. Qiu, A new type of finite volume WENO schemes for hyperbolic conservation laws, J. Sci. Comput., 73 (2017), 1338-1359. [60] J. Zhu and J. Qiu, New finite volume weighted essentially non-oscillatory schemes on triangular meshes, SIAM. J. Sci. Comput., 40(2) (2018), A903-A928. [61] J. Zhu and J. Qiu, A new third order finite volume weighted essentially non-oscillatory scheme on tetrahedral meshes, J. Comput. Phys., 349 (2017), 220-232. [62] J. Zhu and C.-W. Shu, A new type of multi-resolution WENO schemes with increasingly higher order of accuracy, J. Comput. Phys., 375 (2018), 659-683. [63] J. Zhu and C.-W. Shu, A new type of multi-resolution WENO schemes with increasingly higher order of accuracy on triangular meshes, J. Comput. Phys., 392 (2019), 19-33. [64] J. Zhu, X. Zhong, C.-W. Shu and J. Qiu, Runge-Kutta discontinuous Galerkin method using a new type of WENO limiters on unstructured meshes, J. Comput. Phys., 248 (2013), 200-220.

24