Towards optimal high-order compact schemes for simulating compressible flows

Huaibao Zhang; Fan Zhang; Chunguang Xu

doi:10.1016/j.amc.2019.03.001

Towards optimal high-order compact schemes for simulating compressible flows

Zhang, Huaibao;Zhang, Fan;Xu, Chunguang 2018-10-16 00:00:00 Weighted compact nonlinear schemes (WCNS) [Deng and Zhang, JCP 165(2000): 22-44] were developed to improve the performance of the compact high-order non- linear schemes (CNS) by utilizing the weighting technique originally designed for WENO schemes, and excellent shock capturing capability and high resolution are achieved. Various work has been given for further improving the performance of WC- NSs since then. In this work, the ENO-like stencil selection procedure of Targeted ENO schemes [Fu et al. JCP 305(2016):333-359] is introduced for interpolating mid- point variables, targeting compact nonlinear schemes which fully abandon the oscilla- tory stencils crossing discontinuities, and directly apply optimal linear weights when the ﬂow ﬁeld is smooth, such that the optimal numerical resolution is fully recovered in smooth ﬂow ﬁeld. Several canonical numerical cases of scalar equations and the Eu- ler equations of gas dynamics are given to examine the performance of the presented method. Keywords: compact nonlinear schemes; ENO-like stencil-selection; optimal linear weights; gas dynamics 1. Introduction The spatial solution of ﬂow ﬁeld containing strong discontinuities such as shock waves is a challenging topic ever since the shock-capturing schemes invented. High- Corresponding author. Email addresses: zhanghb28@mail.sysu.edu.cn (Huaibao Zhang), zhangfan3@sysu.edu.cn, a04051127@mail.dlut.edu.cn (Fan Zhang) Preprint submitted to Elsevier October 19, 2018 arXiv:1810.07806v1 [physics.comp-ph] 16 Oct 2018 order accurate and high-resolution schemes with discontinuity-capturing ability are especially desired in simulating multi-scale compressible ﬂows. To achieve higher- order accuracy in smooth ﬂow regions, and to suppress spurious oscillations near dis- continuities are two essential problems which especially attract the attention of re- searchers. These two problems, however, somehow contradict with each other, since oscillation usually comes with the use of high-order reconstruction schemes, and high- order schemes owning discontinuity-capturing capability tends to locally degenerate to lower-order schemes for robustness or stability purpose. Great efforts have been made in the development of new high-order discontinuity- capturing schemes to address the dilemma. Among those schemes, the Weighted non- linear schemes, speciﬁcally Weighted Essentially Non-Oscillatory (WENO) schemes [1, 2, 3, 4] and Weighted Compact Nonlinear Schemes (WCNS) [5] are two typical classes of methods, which have achieved great success. They can adaptively tune the numerical dissipation by changing the nonlinear convex combination of candidate sten- cils according to local ﬂow features, such that high order accuracy and non-oscillatory property near discontinuities can be achieved. Compact ﬁnite difference schemes have shown spectral-like resolution [6], which is highly favored in the simulations of multi-scale ﬂow problems. In certain circum- stances, WCNS also has several advantages over the standard ﬁnite-difference WENO: (1) the resolution is slightly higher; (2) various numerical ﬂux schemes can be applied, including Roe’s ﬂux difference splitting (FDS) scheme [7], van Leer’s ﬂux vector split- ting scheme [8], and Liou’s advection upstream splitting method (AUSM) [9]; and (3) WCNS performs well on freestream and vortex preservation properties on wavy grids [10]. Therefore, various researches have been attracted in the improvements and appli- cations of WCNS. The complete WCNS procedure includes three steps [5]: (1) node- to-midpoint weighted nonlinear interpolation of given variables, (2) ﬂux evaluation at the midpoint, and (3) midpoint-to-node high-order differencing of ﬂux function. In the second step, various upwind ﬂux schemes can be applied, and related investigations have been given [11, 12, 13]. In the third step, various candidate midpoint-to-node or midpoint-and-node-to-node difference schemes are also available. Deng and Zhang 2 [14] indicate that the type of midpoint-to-node difference scheme, explicit or implicit, does not affect the resolution of a ﬁfth- or fourth-order WCNS signiﬁcantly, because the weighted nonlinear interpolation, i.e. step (1), dominates the resolution property. In addition, Nonomura et al. [15] demonstrate that the type of midpoint-to-node differ- ence scheme does not signiﬁcantly change the resolution, even for higher-order WC- NSs. Moreover, Nonomura and Fujii [16] proposed an explicit formula which can signiﬁcantly improve the robustness of WCNS and this formula is also more compact because it is featured by using a midpoint-and-node-to-node difference scheme. In gen- eral, aforementioned work suggests that explicit difference scheme is preferred, due to satisfying resolution with lower cost and simplicity for parallelization and vectoriza- tion. Resolutions of WCNSs have been greatly improved by focusing on the improve- ment of its ﬁrst construction step, i.e. node-to-midpoint weighted nonlinear interpola- tion. Novel nonlinear weights achieving optimal order of accuracy [3, 4] are introduced into WCNSs [17, 18], and novel nonlinear compact node-to-midpoint interpolation is also implemented [19]. Recently, a family of high-order targeted ENO schemes has been proposed [20]. Apart from using novel nonlinear weights and incremental-width candidate stencils, the TENO scheme is characterized by using a newly developed ENO-like stencil-selection procedure. These features of TENO schemes bring sig- niﬁcant improvement in spatial resolution. Particularly, the ENO-like stencil-selection procedure is essential to fully recover the background linear schemes in smooth region. Further developments and applications of TENO schemes can be found in [21, 22, 23, 24], and boundary variation diminishing (BVD) method was also used to sharpen the captured shock waves with using TENO scheme [25]. In this work, which is the fur- ther development of compact nonlinear schemes [26], the speciﬁc ENO-like stencil- selection procedure is adopted with an aim that the optimal node-to-midpoint interpo- lation of compact linear scheme is achieved in smooth ﬁelds, including in the region of smooth critical points. This article is organized as follows. In section 2, the WCNS scheme and the pre- sented method, labeled as TCNS similarly to the nomenclature of TENO, are intro- duced in detail. In section 3, detailed numerical analysis including ADR (Approxi- 3 mate dispersion relation) analysis [27] and the solutions of scalar equation and system equations, are introduced to testify the performance of the presented method. Finally, concluding remarks are given in the last section. 2. Numerical methods The governing equations of compressible ﬂows are hyperbolic systems. With- out loss of generality, the theoretical analysis and numerical solutions of the one- dimensional scalar hyperbolic conservation law can be ﬁrst used to examine the perfor- mance of numerical schemes, and then the associated results can be extended to one- or two-dimensional hyperbolic system of equations without substantial difﬁculty. The one-dimensional hyperbolic conservation law can be written as ¶ u ¶ f (u) + = 0; (1) ¶ t ¶ x ¶ f (u) in which the characteristic velocity is and assumed to be positive, without loss ¶ u of generality. Here, the spatial discretization of Eq. (1) is given on an equally spaced one-dimensional mesh, in which the distance between two grid nodes are h, leading to an ordinary differential equation (ODE) system, i.e. du ¶ f = j = f ; i = 1; ; n: (2) x=x i dt ¶ x The numerical solutions of f will be discussed in the following sections, and the tem- poral solutions are given by using the third-order strongly stable Runge-Kutta method [28]. 2.1. The approximation of the ﬂux derivative The ﬁrst-order derivative of ﬂux function, i.e. f , can be approximated by using various numerical schemes. In a WCNS-type method, usually a central (compact) difference scheme is used to calculate the ﬂux derivatives at grid nodes, by utilizing the ﬂux function at midpoints and grid nodes. As aforementioned, because of their efﬁciency, in this work, explicit central dif- ference schemes which require nontridiagonal inversion are used to approximate the 4 derivative. Here, for achieving ﬁfth-order overall accuracy, three midpoint-to-node (MD) or midpoint-and-node-to-node (MND) formulas of sixth-order accuracy are given as a a a 1 2 3 e e e e e e f = ( f f )+ ( f f )+ ( f f ); (3) 1 1 3 3 5 5 i+ i i+ i i+ i 2 2 2 2 2 2 h h h b b b 1 2 3 e e e e e e f = ( f 1 f 1 )+ ( f f )+ ( f f ); (4) i+1 i1 i+2 i2 i+ i 2 2 h h h and c c c 1 2 3 e e e e e e f = ( f 1 f 1 )+ ( f f )+ ( f 3 f 3 ); (5) i+1 i1 i+ i i+ i h 2 2 h h 2 2 which are corresponding to the WCNS [29], Hybrid cell-edge and cell-node Weighted Compact Nonlinear Scheme (HWCNS) [30] and WCNS-midpoint-and-node-to-node difference (WCNS-MND) [16], respectively. Moreover, these three formulas show that they are not compact schemes, but can be taken as the special cases of the general WCNS [29]. The constant coefﬁcients in Eq. (3), (4) and (5), can be found in the corresponding references. It is worth noticed that the required midpoint stencils of HWCNS and WCNS- MND are more compact than that of Eq.(3), because they use node variables in the difference schemes, involving fewer midpoint interpolations. Moreover, Nonomura and Fujii [16] proved that WCNS-MND is more robust but less accurate than the orig- inal WCNS [5]and the original WCNS has higher propensity to blow up compared to WENO schemes when strong discontinuities are captured. WCNS-MND method is thus used in simulations of turbulence problems [31] and multi-species ﬂow problems [32]. In this work, most of the numerical results and discussions are given by using the sixth-order explicit MD scheme (Eq. (3)), and, of course, the other two formulas, Eq. (4) and (5), can be used straightforwardly. 2.2. Node-to-midpoint interpolation: WCNS As shown in the last subsection, midpoint ﬂux terms are unknown and should be evaluated before performing the MD or MND procedure. To achieve this goal, numer- 5 ical upwind ﬂux functions are used. Without loss of generality, a scalar form is used to explain the interpolation procedure. The scalar upwind ﬂux function is h i f = f (u )+ f (u ) ja ˆj u u ; (6) 1 1 1 1 1 i R;i L;i R;i L;i 2 2 2 2 2 where the subscript L and R indicate the variables at the left and right side of x , respectively, and a ˆ is the approximate eigenvalue. The question left to be settled is how to calculate u accurately, while maintaining numerical stability, nonoscillatory L=R;i and sharp discontinuity-capturing properties. For simplicity, we only consider the eval- uation of variables on the left side of x 1 , i.e. u 1 , in the following paragraphs. The i+ L;i+ 2 2 interpolations of the other midpoint variables are performed by using a symmetrical form of u 1 . L;i+ The ﬁfth-order node-to-midpoint reconstruction of WCNS to be used in this work was introduced in [5]. The basic idea of it is to approximate a linear optimal approxi- mation of the midpoint variable u 1 = u + (3u 20u 38u + 60u 5u ); (7) i i2 i1 i i+1 i+2 L;i+ 2 128 by using a nonlinear convex combination of lower-order interpolations. It is obvi- ous that a ﬁve-point full stencil S = fx ; x ; x ; x ; x g is used to calculate i2 i1 i i+1 i+2 i+ the high-order approximation. This optimal ﬁfth-order scheme can be equivalently represented by using three third-order polynomials each constructed on the following three-point substencil S =fx ; x ; x g; k = 1; 2; 3: (8) i+k3 i+k2 i+k1 i+ ;k Each of the third-order polynomials can be expressed in a generic form using the (approximated) nth derivatives (n = 1; 2) Dx (1) (2) u = u (x +Dx) = u + u Dx+ u ; (9) i i i i;k i;k L;i+ ;k where Dx = x x = . Speciﬁcally, the ﬁrst- and second-order derivatives are i+ 2 (1) u = (u 4u + 3u ); i2 i1 i i;1 2h (1) (10) u = (u u ); i+1 i1 i;2 2h (1) u = (3u + 4u u ); i i+1 i+2 i;3 2h 6 and (2) u = (u 2u + u ); i2 i1 i i;1 (2) (11) u = (u 2u + u ); i1 i i+1 i;2 (2) u = (u 2u + u ); i i+1 i+2 i;3 respectively. The linear optimal scheme is then represented as u = d u ; (12) 1 1 å k L;i+ L;i+ ;k 2 2 k=1 where the optimal linear weights are 1 10 5 d = ; d = ; d = : (13) 1 2 3 16 16 16 Obviously, directly using the optimal weights yields a ﬁfth-order accurate scheme. The resulting scheme, however, is inadequate to obtain a satisfying resolution for a steep gradient solution: it causes spurious oscillations. The nonlinear weights of Jiang and Shu [2] can be applied to replace the optimal linear weights in order to suppress non-physical oscillations. The so called JS nonlinear weights are deﬁned by a d k k w = ; a = ; (14) k k 3 2 (b +e) å k k=1 where the small parameter e = 10 is speciﬁed to avoid the denominator becoming zero, and b is the smoothness indicator, which is deﬁned as 2 2 (1) (2) b = hu + h u : (15) i;k i;k The main advantage of using nonlinear weights is that the contribution of oscillatory stencil(s) will be approximately eliminated in the ﬁnal interpolation while preserving those of relative smooth, and thus the spurious numerical oscillation can be suppressed. 2.3. Node-to-midpoint interpolation: ENO-like stencil-selection Instead of merely concentrating on improving the nonlinear weights, the ENO-like stencil-selection procedure [20] is introduced as an essential component of the pre- sented method. Firstly, the nonlinear smoothness measurement yielding strong scale- separation mechanism of the ﬁfth-order scheme is given as g = C + ; k = 1; 2; 3: (16) b +e 7 where t = jb b j is the global smooth indicator that was originally introduced in 5 1 3 [4], and the small threshold is deﬁned as e = 10 . Constant C = 1 is set, and the integer power q = 6 is used. As introduced by Fu et al. [20], larger integer power exponent q and smaller C are preferable for a stronger separation between resolved and non-resolved scales, and the discontinuity-detection capability can be signiﬁcantly enhanced. Rather than using the original nonlinear smoothness measurement, Eq. (16) is fur- ther normalized by c = ; (17) k=1 which is then subject to a sharp cutoff function 0; if c < C ; k T d = (18) 1; otherwise: The underlying idea of this design is that by introducing a value C as the threshold of smoothness, the candidate stencils are attributed as “smooth” or “oscillatory”, such that those stencils of oscillations are abandoned thereby, and only smooth ones are saved with their associated optimal linear weights used in the ﬁnal interpolation. The resulting weight functions are now given by d d (T) k k w = ; (19) d d k k k=1 where d denotes the optimal linear weight, and the above cut-off function d is incor- k k porated into the ﬁnal weighting to decide whether each candidate stencil is taken into account or not. Two promising advantages of this weighting procedure can be readily shown. Firstly, by removing the stencil crossing discontinuity completely, it ensures the numerical robustness of the scheme. Secondly, the background ﬁfth-order linear scheme can be fully recovered in smooth regions, including at smooth critical points. In fact, the stencil-selection procedure has simpliﬁed the nonlinear convex combi- nations of candidate stencils. Instead of using the continuous varying weights lead- ing to inﬁnite possible combinations, the presented method in fact uses only sev- eral candidate combinations, since each candidate stencil has only two possible re- 8 Table 1: The coefﬁcients of the equivalent single polynomial spatial reconstructions. if d = u ˆ S a a a a a 1;2;3 m;i2 m;i1 m;i m;i+1 m;i+2 1 m L;i+ ;m 1,1,1 u ˆ S 3/128 -5/32 45/64 15/32 -5/128 1 0 L;i+ ;0 0,1,1 u ˆ S 0 -1/12 5/8 1/2 -1/24 1 1 L;i+ ;1 1,1,0 u ˆ S 3/88 -5/22 75/88 15/44 0 1 2 L;i+ ;2 0,0,1 u ˆ S 0 0 3/8 3/4 -1/8 1 3 L;i+ ;3 0,1,0 u ˆ S 0 -1/8 3/4 3/8 0 1 4 L;i+ ;4 1,0,0 u ˆ S 3/8 -5/4 15/8 0 0 1 5 L;i+ ;5 1,0,1 u ˆ S 1/16 -5/24 5/8 5/8 -5/48 1 6 L;i+ ;6 Figure 1: Schematic of the equivalent candidate stencils of the ﬁfth-order nonlinear interpolation. The circle indicates that the stencil is not really continuous. sults of d d , i.e. d and 0. Therefore, similar to [24], each of the possible com- k k k binations of the nonlinear interpolation of u ˆ 1 can be represented as u ˆ = L;i+ L;i+ ;m a u + a u + a u + a u + a u , where the coefﬁcients are m;i2 i2 m;i1 i1 m;i i m;i+1 i+1 m;i+2 i+2 given in Table 1. The corresponding equivalent stencil of each convex combination is given in Fig.1, by which the meaning of ”stencil-selection” can be interpreted more clear. In addition, using the stencil-selection procedure results into several potential choices of high-order polynomials, and each can be optimized independently such as in [24]. Parameter C , serving as the global reference of smooth indicators, is an effective and a direct mean to control the spectral properties of the scheme for compressible turbu- lence simulation in which close embedded shocklets need to be captured [24]. Within 9 5 the scope of this work, however, the parameter is simply set as C = 10 for all the simulations without further investigation on the inﬂuence of its choice, but the advan- tages of importing the ENO-like stencil-selection procedure into compact nonlinear schemes are clearly showed. 3. Numerical results A variety of canonical problems are simulated to assess ﬁfth-order WCNS-JS, WCNS-Z, and the proposed scheme TCNS. One-dimensional linear advection equa- tion, one-dimensional inviscid Burgers equation and Euler equations of gas dynamics are used as model equations. The ideal-gas equation of state is given by p = (g 1)r e with g = 1:4 to close Euler equations. Node-to-midpoint interpolation is performed on characteristic variables to alleviate spurious oscillations [5]. Van leer scheme [8] is used for the computation of ﬂuxes. The proposed ﬁfth-order TCNS scheme is ﬁrst assessed through comparisons with the classical WCNS-JS, WCNS-Z on spectral prop- erties and ofﬁcial order of accuracy. Its rate of convergence is then veriﬁed by solv- ing one-dimensional linear advection equation, and one-dimensional inviscid Burgers equation. Performances of the above ﬁfth-order schemes in shock-capturing and wave resolutions are compared by solving selected test-problems both in one and two dimen- sions. The CFL= 0:6 is used as default for all numerical schemes and testing cases. 3.1. Approximate dispersion relation analysis Following the ADR analysis introduced by Pirozzoli [33] and Tu et al. [34], the spectral properties of different ﬁfth-order schemes are compared and shown in Fig. 2. The result of the presented method agrees very well with the background linear scheme for the low range of wave-numbers with a recovered wave-number reaching up to 1.76. A signiﬁcant improvement can also be found when cosmpared against the WCNS-JS in both dispersion and dissipation properties. While slight inferior result of dissipation property is given by TCNS when compared with WCNS-Z for the wave-number in intermediate range, TCNS performs better than WCNS-Z scheme for rest range of wave-numbers. Moreover, the performance of TCNS will be signiﬁcantly improved if 10 small C is used, and the adaptive dissipation control in [24] will further improve its overall performance. (a) (b) Figure 2: Approximated dispersion and dissipation properties of different ﬁfth-order schemes: dispersion (left) and dissipation (right). 3.2. Linear advection problem The one-dimensional Gaussian pulse advection problem [35] is used to assess the numerical order of accuracy of the proposed scheme. This problem is modeled by means of the linear advection equation, given by ¶ u ¶ u + = 0; x2 [0; 1]; (20) ¶ t ¶ x with periodic boundary-conditions and the initial conditions given by u(t; x = 0) = u(t; x = 1); (21) 300(xx ) u(0; x) = e ; x = 0:5: Time integration is performed up to t = 1, which corresponds to one period of the single wave propagation in time. A set of evenly distributed grids are progressively reﬁned by a factor of 2 from the most coarse grid with N = 51. Numerical simulation on each of the grid is conducted by using a time step size that even further reducing its value will not change the evaluated error of the numerical solution. 11 Table 2 and Fig. 3 illustrate the numerical error of those ﬁfth-order schemes used, where the result of TCNS scheme coincides with that of the corresponding linear scheme, indicating that it recovers the optimal order of convergence. A slight devi- ation can be found for WCNS-Z over the linear scheme for the coarse grids. WCNS-JS also shows approximate ﬁfth-order accuracy, but its resolution is signiﬁcantly lower than WCNS-Z and TCNS. In general, using the ENO-like stencil-selection procedure recovers the optimal linear scheme in this smooth ﬁeld. Table 2: L -error and convergence rate for different ﬁfth-order schemes used by solving the linear advection equation at t = 1 . Linear WCNS-JS WCNS-Z TCNS Error Order Error Order Error Order Error Order 51 5.22E-02 * 1.07E-01 * 5.67E-02 * 5.20E-02 * 101 3.30E-03 3.98 1.04E-02 3.37 3.57E-03 3.99 3.30E-03 3.98 201 1.16E-04 4.83 4.63E-04 4.49 1.20E-04 4.90 1.16E-04 4.83 401 3.69E-06 4.97 1.84E-05 4.66 3.72E-06 5.01 3.69E-06 4.97 801 1.16E-07 4.99 6.36E-07 4.85 1.16E-07 5.00 1.16E-07 4.99 1601 3.64E-09 4.99 2.02E-08 4.98 3.64E-09 4.99 3.64E-09 4.99 (a) L error (b) L error 2 ¥ Figure 3: Convergence rate of the L - and L -error for different ﬁfth-order schemes used by solving the 2 ¥ linear advection equation at t = 1. 12 3.3. Inviscid Burgers equation The one-dimensional inviscid Burgers equation [36, 37] is used to assess the actual order of accuracy of the proposed scheme when it is applied to a non-linear scalar equation. The governing equation is in the form of ¶ u ¶ 1 + u = 0; x2 [0; 2]; (22) ¶ t ¶ x 2 with periodic boundary-conditions and the initial conditions given by u(t; x = 0) = u(t; x = 1); (23) u(0; x) = + sin(p x): Exact solutions are computed by solving the derived general characteristic relation in reference [36]. The solutions are smooth for 0 t < 1=p , and a discontinuity devel- ops and starts to interact with the expansion wave if t 1=p . The results at t = 0:2 and t = 0:7 are both given to show the continuous and discontinuous distributions. Solu- tions obtained by the above ﬁfth-order schemes are compared against the exact solution in Fig. 4, in which, good agreement is presented both for the smooth and the discon- tinuous solution. Moreover, at t = 0:2, L -error and convergence rate for each scheme are presented in Table 3. All of the schemes can achieve ﬁfth-order of accuracy as the grids are reﬁned. Perfect agreement can be found between the TCNS and the back- ground ﬁfth-order linear scheme. WCNS-Z scheme shows minor deviation compared with TCNS scheme only when the grid is relatively coarse. WCNS-JS scheme is less accurate than the other ﬁfth-order schemes. Again, TCNS scheme recovers the optimal linear scheme if the solution is smooth. Moreover, it should be noticed that in the discontinuous solution, minor overshoot has been found. This problem can be avoided by using the MND method, as shown in Fig.4(b). At the meantime, MND method reduces the resolution of the result, as already introduced by Nonomura and Fujii [16]. Whereas, in this work, only the node- to-midpoint interpolation methods are discussed, and thus the detail of MND method is not further investigated. As mentioned above, the presented method can be applied with using all the midpoint-to-node or midpoint-and-node-to-node methods. 13 (a) t = 0:2 (b) t = 0:7 Figure 4: One-dimensional inviscid Burgers equation: numerical solutions and exact solution at t = 0:2 (left) and t = 0:7 (right). Grid resolution 161. Table 3: L -error and convergence rate for different ﬁfth-order schemes used by solving the 1-D inviscid Burgers equation at t = 0.2. Linear WCNS-JS WCNS-Z TCNS Error Order Error Order Error Order Error Order 41 9.96E-04 * 1.41E-03 * 1.01E-03 * 9.96E-04 * 81 9.04E-05 3.46 1.34E-04 3.40 9.07E-05 3.48 9.04E-05 3.46 161 3.21E-06 4.81 5.06E-06 4.73 3.22E-06 4.82 3.21E-06 4.81 321 1.00E-07 5.00 1.63E-07 4.96 1.00E-07 5.00 1.00E-07 5.00 641 2.82E-09 5.15 4.82E-09 5.07 2.82E-09 5.15 2.82E-09 5.15 3.4. Sod and Lax shock tube problem Riemann initial value problems of Sod [38] and Lax [39] are used to evaluate the shock-capturing capability of the proposed scheme. The initial conditions for the Sod problem is (1; 0; 1) x2 [0; 0:5]; (r; u; p) = (24) (0:125; 0; 0:1) x2 (0:5; 1]: and the ﬁnal time is t = 2 by solving the problem on an evenly distributed grid of N = 101 points. 14 (a) (b) Figure 5: Sod problem: numerical solutions and the exact solution at t = 0:2. Overall, good agreement of numerical solutions with the exact solution is achieved for all of the schemes shown in Fig. 5(a). Steep gradient solutions such as shock wave and contact discontinuity are captured without oscillations. In Fig. 5(b), WCNS-Z and TCNS schemes yield a relative sharp proﬁle compared to WCNS-JS in the vicinity of contact discontinuity, and TCNS is found with a better performance of resolving the contact discontinuity when compared with WCNS-Z. The initial conditions for the Lax shock-tube problem is (0:445; 0:698; 3:528) x2 [0; 0:5]; (r; u; p) = (25) (0:5; 0; 0:571) x2 (0:5; 1]: This case is run on a grid of N = 101 points with uniform distribution, and the results at t = 0:14 are given. Density and velocity distributions are compared in Fig.6(a), and Fig. 6(b) respec- tively. WCNS-Z overshoots the velocity slightly at the tail of the expansion fan, but TCNS shows accurate and oscillation-free result. WCNS-JS smears the shock wave signiﬁcantly, indicating extra numerical dissipation. 15 (a) Density (b) Velocity Figure 6: Lax problem: numerical solutions and the exact solution at t = 0:14. 3.5. Shock-density wave interaction The shock-density wave interaction problem [40] is characterized by a right mov- ing Mach 3 shock interacting with sine waves in density ﬁeld. The multi-scale wave structure is evolved after the shock wave interacts with the oscillating density wave, and both the shock-capturing and wave-resolution capabilities are evaluated thereafter. The problem is initialized by (3:857; 2:629; 10:333); x2 [0; 1]; (r; u; p) = (26) (1+ 0:2sin(5x); 0; 1); x2 (1; 10]: This case is run on a grid of N = 201 points which are uniformly distributed and the ﬁnal time is t = 1:8. Numerical solution of WCNS-JS on a grid of N = 2001 is used as the reference ”exact” solution. As shown in Fig. 7, TCNS produces considerably better resolved density waves be- hind the shock wave compared with WCNS-JS and WCNS-Z. Particularly, the result of WCNS-JS indicates strong numerical dissipation, since the small-scale wave structure is signiﬁcantly smeared. Moreover, TCNS does not produce spurious oscillations in the ﬂow ﬁeld. 16 (a) (b) Figure 7: Shock-density wave interaction problem: numerical solutions and the exact solution at t = 1:8. 3.6. Two-dimensional Riemann problem In this section, the conﬁgurations 3 and 6 out of 19 2-D Riemann problems used by Lax and Liu [41] are speciﬁcally taken to evaluate performances of the proposed numerical schemes. 3.6.1. The conﬁguration 3 The initial conditions for the conﬁguration 3 are given by 1 1 (1:5; 0; 0; 1:5) x2 ; 1 and y2 ; 1 ; > 2 2 1 1 (0:5323; 1:206; 0; 0:3) x2 0; and y2 ; 1 ; 2 2 (r; u; v; p) = (27) > 1 1 >(0:138; 1:206; 1:206; 0:029) x2 0; and y2 0; ; 2 2 : 1 1 (0:5323; 0; 1:206; 0:3) x2 ; 1 and y2 0; : 2 2 Boundary conditions are given by ¶u(t; x; y) = 0; x = 0; 1; 8t; y; ¶ x (28) ¶u(t; x; y) = 0; y = 0; 1; 8t; x: ¶ y Figure 8 presents the results obtained by using WCNS-JS, WCNS-Z and TCNS schemes, at t = 0:3. The density proﬁles are best predicted by TCNS. WCNS-Z yields 17 a solution very close to that of TCNS, but TCNS captures more vortexes along the contact lines. WCNS-JS smears a few small scales that can be resolved by TCNS and WCNS-Z due to relative large dissipations. (a) WCNS-JS (b) WCNS-Z (c) TCNS Figure 8: Conﬁgurations 3 of 2-D Riemann problems in [41]: 30 density contours ranging from 0.1 to 1.8 at t = 0:23 obtained on a grid of 1024 2014. 18 3.6.2. The conﬁguration 6 The initial conditions for the conﬁguration 6 are given by 1 1 (1; 0:75;0:5; 1) x2 ; 1 and y2 ; 1 ; > 2 2 1 1 (2; 0:75; 0:5; 1) x2 0; and y2 ; 1 ; 2 2 (r; u; v; p) = (29) 1 1 > (1;0:75; 0:5; 1) x2 0; and y2 0; ; 2 2 : 1 1 (3;0:75;0:5; 1) x2 ; 1 and y2 0; : 2 2 Boundary conditions are the same as in the last case. Results obtained by using WCNS-JS, WCNS-Z and TCNS schemes at t = 0:3 are shown in Fig. 9. The results of WCNS-JS and WCNS-Z are similar. Whereas, using TCNS obtains abundant small scale ﬂow structures along the contact lines, indicating lower numerical dissipation. 19 (a) WCNS-JS (b) WCNS-Z (c) TCNS Figure 9: Conﬁgurations 6 of 2-D Riemann problems in [41]: 40 density contours ranging from 0.1 to 2.9 at t = 0:3 obtained on a grid of 1024 2014. 3.7. Shock vortex interaction This problem was modeled in two dimensions by Jiang and Shu [2]. It involves a vortex perturbing a stationary shock. The computational domain is [0; 2] [0; 1]. This ﬁeld is discretized by 251 101 grid. Initially, a stationary shock is located at x = 0:5 normal to the x direction. The left state of this shock is speciﬁed as (r; u; v; p) = (1; 1:1 g; 0; 1). A small vortex centered at (0.25, 0.5) is superposed to the ﬂowﬁled on the left hand side of the normal shock. The superposition is performed through 20 perturbations on velocity (u; v), temperature T , given by T = p=r , and entropy S, deﬁned as S = ln( p=r ) of the mean ﬂow. Speciﬁcally, the perturbation variables are a(1t ) u ˜ = et e sinq a(1t ) u ˜ =et e cosq (30) 2 2a(1g ) (g1)e e > T = 4ag S = 0 2 2 where t = r=r , and r = (x x ) +(y y ) , r = 0:05, e = 0:3, a = 0:204 are taken c c c c from the reference [2]. The results at t = 0:6 are presented in Fig.10. (a) WCNS-JS (b) WCNS-Z (c) TCNS Figure 10: Shock vortex interaction problem: 90 density contours ranging from 1.19 to 1.37. 21 The ﬂow contours of different schemes are similar. The density distribution along the central line of the ﬂow ﬁeld are further shown in Fig.11, where the result obtained by using WCNS-Z on a reﬁned mesh of 1001 401 is used as the “exact” solution. It can be found that the captured vortex of TCNS is more accurate, comparing with the other two results of WCNS-JS and WCNS-Z. Figure 11: Shock vortex interaction problem: density distribution along y = 0:5. 3.8. Rayleigh-Taylor instability Rayleigh-Taylor instability problem which contains both discontinuities and com- plex ﬂow structures [42], is also used to examine the numerical dissipation of the pre- sented method. The initial conditions are given by (2; 0;0:025a cos(8p x); 1+ 2y) x2 [0; 0:25] and y2 [0; 0:5); (r; u; v; p) = (1; 0;0:025a cos(8p x); 1+ 3=2) x2 [0; 0:25] and y2 [0:5; 1]; (31) where a is the speed of sound, given by a = g and a different g = is used for r 3 this speciﬁc case. Reﬂecting boundary conditions are imposed at the left and right side of the domain, and constant boundary conditions are given for the top and the bottom 22 sides, in details (1; 0; 0; 2:5) y = 1; 8t; x; (r; u; v; p) = (32) (2; 0; 0; 1) y = 0; 8t; x: Two source terms r , and r v are added to the right hand side of the third and the fourth equation, respectively. Two sets of grids are used, i.e. 128 512 and 256 1024. Density proﬁles at ﬁnal time t = 1:95 calculated by WCNS-JS, WCNS-Z and TCNS are shown in Fig. 12. 23 (a) Grid resolution 128 512 (b) Grid resolution 256 1024 Figure 12: Rayleigh-Taylor instability problem: 30 density contours ranging from 0.9 to 2.2. 24 TCNS has resolved much more abundant vortical structures than WCNS-JS and WCNS-Z, suggesting that TCNS is signiﬁcantly less dissipative. Furthermore, the low-dissipation property of TCNS induces symmetry breaking of the ﬂow ﬁeld. In par- ticular, the small structures resolved by using TCNS with grid resolution of 128 512 are comparable to that of using WCNS-JS with grid resolution of 256 1024. Even the result of WCNS-Z only shows little advantage (if any) while doubling the grid resolution. 4. Conclusions As a subsequent work of compact nonlinear scheme and weighed compact non- linear scheme, we have introduced a novel compact nonlinear scheme in this article, by applying the ENO-like stencil-selection procedure. The method is targeted optimal linear interpolation in the node-to-midpoint interpolation step of the compact nonlin- ear scheme, and thus the method is named as targeted CNS (TCNS). Numerical results of one-dimensional scalar equations and one- or two-dimensional Euler equations are given to examine the performance of the method, involving strong discontinuities and broadband ﬂuctuations. The presented method uses a novel smoothness measurement which possesses strong scale separation mechanism, and the ENO-like stencil-selection procedure is capable to target optimal linear weights while maintaining excellent shock-capturing capability. ADR analysis shows that TCNS recovers the background linear scheme up to higher wave number, even using a relative larger C . Examining the numerical results, it can be found that TCNS has recovered optimal linear scheme in smooth ﬁeld, where the optimal linear weights d is directly applied in the node-to-midpoint interpolation. Whereas, WCNSs, including WCNS-JS and WCNS-Z, approximate the optimal linear weights but can not perfectly recover them. Therefore, the presented TCNS is capa- ble to capture much more abundant wave structures in the simulations. Especially in the Rayleigh-Taylor instability problem, TCNS achieves similar or even better result comparing with WCNS-JS on a coarser grid. 25 Acknowledgments Financial support for this work was provided through grant number 20187413071020008. References [1] X. D. Liu, S. Osher, T. Chan, Weighted essentially non-oscillatory schemes, Jour- nal of Computational Physics 115 (1) (1994) 200–212. [2] G.-S. Jiang, C.-W. Shu, Efﬁcient implementation of weighted ENO schemes, Journal of Computational Physics 126 (1) (1996) 202–228. [3] A. K. Henrick, T. D. Aslam, J. M. Powers, Mapped weighted essentially non- oscillatory schemes: Achieving optimal order near critical points, Journal of Computational Physics 207 (2) (2005) 542–567. [4] R. Borges, M. Carmona, B. Costa, W. S. Don, An improved weighted essentially non-oscillatory scheme for hyperbolic conservation laws, Journal of Computa- tional Physics 227 (6) (2008) 3191–3211. doi:10.1016/j.jcp.2007.11.038. [5] X. Deng, H. Zhang, Developing high-order weighted compact nonlinear schemes, Journal of Computational Physics 165 (1) (2000) 22 – 44. doi:10.1006/jcph. 2000.6594. [6] S. K. Lele, Compact ﬁnite difference schemes with spectral-like resolution, Journal of Computational Physics 103 (1) (1992) 16 – 42. doi:10.1016/ 0021-9991(92)90324-R. [7] P. L. Roe, Approximate Riemann solvers, parameter vectors, and difference schemes, Journal of Computational Physics 43 (2) (1981) 357–372. [8] B. van Leer, Flux-vector splitting for the Euler equations, in: E. Krause (Ed.), Eighth International Conference on Numerical Methods in Fluid Dynamics, Vol. 170 of Lecture Notes in Physics, Springer Berlin Heidelberg, Aachen, Germany, 1982, pp. 507–512. 26 [9] M. S. Liou, On a new class of ﬂux splittings, in: M. Napolitano, F. Sabetta (Eds.), Thirteenth International Conference on Numerical Methods in Fluid Dynamics, Vol. 414 of Lecture Notes in Physics, Springer Berlin Heidelberg, Berlin, Heidel- berg, 1993, pp. 115–119. [10] T. Nonomura, N. Iizuka, K. Fujii, Freestream and vortex preservation properties of high-order weno and wcns on curvilinear grids, Computers & Fluids 39 (2) (2010) 197 – 214. doi:10.1016/j.compfluid.2009.08.005. [11] D. Wang, X. Deng, G. Wang, Y. Dong, Developing a hybrid ﬂux function suitable for hypersonic ﬂow simulation with high-order methods, Int. J. Numer. Meth. Fluids 81 (5) (2016) 309–327. [12] G.-H. TU, J.-Q. CHEN, M.-L. MAO, X.-H. ZHAO, H.-Y. LIU, On the splitting methods of inviscid ﬂuxes for implementing high-order weighted compact non- linear schemes, Applied Mathematics and Mechanics 37 (12) (2016) 1324–1344. doi:10.21656/1000-0887.370518. [13] H. Zhu, Z. Yan, H. Liu, M. Mao, Properties of osher ﬂux with entropy ﬁx in high- order WCNS, Acta Aeronautica et Astronautica Sinica 38 (5) (2017) 120543. [14] X. Deng, X. Liu, M. Mao, H. Zhang, Investigation on weighted compact ﬁfth- order nonlinear scheme and applications to complex ﬂow, in: 17th AIAA Com- putational Fluid Dynamics Conference, no. 0 in Fluid Dynamics and Co-located Conferences, AIAA Paper 2005-5246, 2005. doi:10.2514/6.2005-5246. [15] T. Nonomura, K. Fujii, Effects of difference scheme type in high-order weighted compact nonlinear schemes, Journal of Computational Physics 228 (10) (2009) 3533 – 3539. doi:10.1016/j.jcp.2009.02.018. [16] T. Nonomura, K. Fujii, Robust explicit formulation of weighted compact non- linear scheme, Computers & Fluids 85 (2013) 8 – 18. doi:10.1016/j. compfluid.2012.09.001. 27 [17] Z. Yan, H. Liu, M. Mao, H. Zhu, X. Deng, New nonlinear weights for improving accuracy and resolution of weighted compact nonlinear scheme, Computers & Fluids 127 (2016) 226 – 240. doi:10.1016/j.compfluid.2016.01.005. [18] Z.-G. Yan, H. Liu, Y. Ma, M. Mao, X. Deng, Further improvement of weighted compact nonlinear scheme using compact nonlinear interpolation, Computers & Fluids 156 (2017) 135 – 145. doi:10.1016/j.compfluid.2017.06.028. [19] Y. Ma, H. Liu, Z. Yan, M. Mao, X. Deng, Investigation of compact interpolation based on hwcns schemes, Chinese Journal of Computational Mechanics 32 (3) (2015) 388–393. [20] L. Fu, X. Y. Hu, N. A. Adams, A family of high-order targeted ENO schemes for compressible-ﬂuid simulations, Journal of Computational Physics 305 (2016) 333–359. [21] L. Fu, X. Y. Hu, N. A. Adams, A targeted eno scheme as implicit model for tur- bulent and genuine subgrid scales, Journal of Computational Physics 349 (2017) 97–121. [22] O. Haimovich, S. H. Frankel, Numerical simulations of compressible multicom- ponent and multiphase ﬂow using a high-order targeted ENO (TENO) ﬁnite- volume method, Computers & Fluids 146 (2017) 105–116. doi:10.1016/j. compfluid.2017.01.012. [23] L. Fu, X. Hu, N. A. Adams, A targeted eno scheme as implicit model for turbulent and genuine subgrid scales, Communications in Computational Physics. [24] L. Fu, X. Y. Hu, N. A. Adams, A new class of adaptive high-order tar- geted ENO schemes for hyperbolic conservation laws, Journal of Computational Physicsdoi:10.1016/j.jcp.2018.07.043. [25] Z. Sun, S. Inaba, F. Xiao, Boundary variation diminishing (bvd) reconstruction: A new approach to improve godunov schemes, Journal of Computational Physics 322 (2016) 309 – 325. doi:10.1016/j.jcp.2016.06.051. 28 [26] X. Deng, H. Maekawa, Compact high-order accurate nonlinear schemes, Journal of Computational Physics 130 (1) (1997) 77 – 91. doi:10.1006/jcph.1996. [27] S. Pirozzoli, On the spectral properties of shock-capturing schemes, Journal of Computational Physics 219 (2) (2006) 489 – 497. doi:10.1016/j.jcp.2006. 07.009. [28] S. Gottlieb, C. Shu, E. Tadmor, Strong stability-preserving high- order time discretization methods, SIAM Review 43 (1) (2001) 89– 112. arXiv:https://doi.org/10.1137/S003614450036757X, doi:10.1137/S003614450036757X. URL https://doi.org/10.1137/S003614450036757X [29] X. Deng, High-order accurate dissipative weighted compact nonlinear schemes, Science in China Series A: Mathematics 45 (3) (2002) 356. doi:10.1360/ 02ys9037. [30] X. Deng, New high-order hybrid cell-edge and cell-node weighted compact non- linear schemes, in: 20th AIAA Computational Fluid Dynamics Conference, AIAA Paper 2011-3857, 2011. doi:10.2514/6.2011-3857. [31] G. Zhao, M. Sun, S. Xie, H. Wang, Numerical dissipation control in an adaptive wcns with a new smoothness indicator, Applied Mathematics and Computation 330 (2018) 239 – 253. doi:10.1016/j.amc.2018.01.019. [32] M. L. Wong, S. K. Lele, High-order localized dissipation weighted compact non- linear scheme for shock- and interface-capturing in compressible ﬂows, Journal of Computational Physics 339 (2017) 179 – 209. doi:10.1016/j.jcp.2017. 03.008. [33] S. Pirozzoli, On the spectral properties of shock-capturing schemes, Journal of Computational Physics 219 (2006) 489–497. doi:10.1016/j.jcp.2006.07. 29 [34] G.-H. Tu, X.-G. Deng, M.-L. Mao, Spectral property comparison of ﬁfth-order nonlinear wcns and weno difference schemes, Acta Aerodynamica Sinica 30 (6) (2013) 709–712. [35] N. K. Yamaleev, M. H. Carpenter, A systematic methodology for constructing high-order energy stable weno schemes, Journal of Computational Physics 228 (11) (2009) 4248 – 4272. doi:https://doi.org/10.1016/j.jcp. 2009.03.002. URL http://www.sciencedirect.com/science/article/pii/ S0021999109001132 [36] A. Harten, B. Engquist, S. Osher, S. R. Chakravarthy, Uniformly high order ac- curate essentially non-oscillatory schemes, iii, Journal of Computational Physics 131 (1) (1997) 3 – 47. doi:https://doi.org/10.1006/jcph.1996.5632. URL http://www.sciencedirect.com/science/article/pii/ S0021999196956326 [37] G. Gerolymos, D. Snchal, I. Vallet, Very-high-order weno schemes, Journal of Computational Physics 228 (23) (2009) 8481 – 8524. doi:https://doi.org/10.1016/j.jcp.2009.07.039. URL http://www.sciencedirect.com/science/article/pii/ S0021999109003908 [38] G. A. Sod, A survey of several ﬁnite difference methods for systems of nonlinear hyperbolic conservation laws, Journal of Computational Physics 27 (1) (1978) 1 – 31. doi:https://doi.org/10.1016/0021-9991(78)90023-2. URL http://www.sciencedirect.com/science/article/pii/ [39] P. D. Lax, Weak solutions of nonlinear hyperbolic equations and their numerical computation, Communications on Pure and Applied Mathematics 7 (1) 159–193. arXiv:https://onlinelibrary.wiley.com/doi/pdf/10.1002/cpa. 3160070112, doi:10.1002/cpa.3160070112. 30 URL https://onlinelibrary.wiley.com/doi/abs/10.1002/cpa. [40] C.-W. Shu, S. Osher, Efﬁcient implementation of essentially non-oscillatory shock-capturing schemes, ii, Journal of Computational Physics 83 (1) (1989) 32 – 78. doi:https://doi.org/10.1016/0021-9991(89)90222-2. URL http://www.sciencedirect.com/science/article/pii/ [41] P. Lax, X. Liu, Solution of two-dimensional riemann problems of gas dynam- ics by positive schemes, SIAM Journal on Scientiﬁc Computing 19 (2) (1998) 319–340. arXiv:https://doi.org/10.1137/S1064827595291819, doi: 10.1137/S1064827595291819. URL https://doi.org/10.1137/S1064827595291819 [42] Z. Xu, C.-W. Shu, Anti-diffusive ﬂux corrections for high order ﬁnite difference weno schemes, Journal of Computational Physics 205 (2) (2005) 458 – 485. doi:https://doi.org/10.1016/j.jcp.2004.11.014. URL http://www.sciencedirect.com/science/article/pii/ S0021999104004760 http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Physics arXiv (Cornell University) http://www.deepdyve.com/lp/arxiv-cornell-university/towards-optimal-high-order-compact-schemes-for-simulating-compressible-95OmiHl8RW

Loading next page...

References

References for this paper are not available at this time. We will be adding them shortly, thank you for your patience.

ISSN: 0096-3003
eISSN: ARCH-3341
DOI: 10.1016/j.amc.2019.03.001
Publisher site: See Article on Publisher Site

Abstract

Weighted compact nonlinear schemes (WCNS) [Deng and Zhang, JCP 165(2000): 22-44] were developed to improve the performance of the compact high-order non- linear schemes (CNS) by utilizing the weighting technique originally designed for WENO schemes, and excellent shock capturing capability and high resolution are achieved. Various work has been given for further improving the performance of WC- NSs since then. In this work, the ENO-like stencil selection procedure of Targeted ENO schemes [Fu et al. JCP 305(2016):333-359] is introduced for interpolating mid- point variables, targeting compact nonlinear schemes which fully abandon the oscilla- tory stencils crossing discontinuities, and directly apply optimal linear weights when the ﬂow ﬁeld is smooth, such that the optimal numerical resolution is fully recovered in smooth ﬂow ﬁeld. Several canonical numerical cases of scalar equations and the Eu- ler equations of gas dynamics are given to examine the performance of the presented method. Keywords: compact nonlinear schemes; ENO-like stencil-selection; optimal linear weights; gas dynamics 1. Introduction The spatial solution of ﬂow ﬁeld containing strong discontinuities such as shock waves is a challenging topic ever since the shock-capturing schemes invented. High- Corresponding author. Email addresses: zhanghb28@mail.sysu.edu.cn (Huaibao Zhang), zhangfan3@sysu.edu.cn, a04051127@mail.dlut.edu.cn (Fan Zhang) Preprint submitted to Elsevier October 19, 2018 arXiv:1810.07806v1 [physics.comp-ph] 16 Oct 2018 order accurate and high-resolution schemes with discontinuity-capturing ability are especially desired in simulating multi-scale compressible ﬂows. To achieve higher- order accuracy in smooth ﬂow regions, and to suppress spurious oscillations near dis- continuities are two essential problems which especially attract the attention of re- searchers. These two problems, however, somehow contradict with each other, since oscillation usually comes with the use of high-order reconstruction schemes, and high- order schemes owning discontinuity-capturing capability tends to locally degenerate to lower-order schemes for robustness or stability purpose. Great efforts have been made in the development of new high-order discontinuity- capturing schemes to address the dilemma. Among those schemes, the Weighted non- linear schemes, speciﬁcally Weighted Essentially Non-Oscillatory (WENO) schemes [1, 2, 3, 4] and Weighted Compact Nonlinear Schemes (WCNS) [5] are two typical classes of methods, which have achieved great success. They can adaptively tune the numerical dissipation by changing the nonlinear convex combination of candidate sten- cils according to local ﬂow features, such that high order accuracy and non-oscillatory property near discontinuities can be achieved. Compact ﬁnite difference schemes have shown spectral-like resolution [6], which is highly favored in the simulations of multi-scale ﬂow problems. In certain circum- stances, WCNS also has several advantages over the standard ﬁnite-difference WENO: (1) the resolution is slightly higher; (2) various numerical ﬂux schemes can be applied, including Roe’s ﬂux difference splitting (FDS) scheme [7], van Leer’s ﬂux vector split- ting scheme [8], and Liou’s advection upstream splitting method (AUSM) [9]; and (3) WCNS performs well on freestream and vortex preservation properties on wavy grids [10]. Therefore, various researches have been attracted in the improvements and appli- cations of WCNS. The complete WCNS procedure includes three steps [5]: (1) node- to-midpoint weighted nonlinear interpolation of given variables, (2) ﬂux evaluation at the midpoint, and (3) midpoint-to-node high-order differencing of ﬂux function. In the second step, various upwind ﬂux schemes can be applied, and related investigations have been given [11, 12, 13]. In the third step, various candidate midpoint-to-node or midpoint-and-node-to-node difference schemes are also available. Deng and Zhang 2 [14] indicate that the type of midpoint-to-node difference scheme, explicit or implicit, does not affect the resolution of a ﬁfth- or fourth-order WCNS signiﬁcantly, because the weighted nonlinear interpolation, i.e. step (1), dominates the resolution property. In addition, Nonomura et al. [15] demonstrate that the type of midpoint-to-node differ- ence scheme does not signiﬁcantly change the resolution, even for higher-order WC- NSs. Moreover, Nonomura and Fujii [16] proposed an explicit formula which can signiﬁcantly improve the robustness of WCNS and this formula is also more compact because it is featured by using a midpoint-and-node-to-node difference scheme. In gen- eral, aforementioned work suggests that explicit difference scheme is preferred, due to satisfying resolution with lower cost and simplicity for parallelization and vectoriza- tion. Resolutions of WCNSs have been greatly improved by focusing on the improve- ment of its ﬁrst construction step, i.e. node-to-midpoint weighted nonlinear interpola- tion. Novel nonlinear weights achieving optimal order of accuracy [3, 4] are introduced into WCNSs [17, 18], and novel nonlinear compact node-to-midpoint interpolation is also implemented [19]. Recently, a family of high-order targeted ENO schemes has been proposed [20]. Apart from using novel nonlinear weights and incremental-width candidate stencils, the TENO scheme is characterized by using a newly developed ENO-like stencil-selection procedure. These features of TENO schemes bring sig- niﬁcant improvement in spatial resolution. Particularly, the ENO-like stencil-selection procedure is essential to fully recover the background linear schemes in smooth region. Further developments and applications of TENO schemes can be found in [21, 22, 23, 24], and boundary variation diminishing (BVD) method was also used to sharpen the captured shock waves with using TENO scheme [25]. In this work, which is the fur- ther development of compact nonlinear schemes [26], the speciﬁc ENO-like stencil- selection procedure is adopted with an aim that the optimal node-to-midpoint interpo- lation of compact linear scheme is achieved in smooth ﬁelds, including in the region of smooth critical points. This article is organized as follows. In section 2, the WCNS scheme and the pre- sented method, labeled as TCNS similarly to the nomenclature of TENO, are intro- duced in detail. In section 3, detailed numerical analysis including ADR (Approxi- 3 mate dispersion relation) analysis [27] and the solutions of scalar equation and system equations, are introduced to testify the performance of the presented method. Finally, concluding remarks are given in the last section. 2. Numerical methods The governing equations of compressible ﬂows are hyperbolic systems. With- out loss of generality, the theoretical analysis and numerical solutions of the one- dimensional scalar hyperbolic conservation law can be ﬁrst used to examine the perfor- mance of numerical schemes, and then the associated results can be extended to one- or two-dimensional hyperbolic system of equations without substantial difﬁculty. The one-dimensional hyperbolic conservation law can be written as ¶ u ¶ f (u) + = 0; (1) ¶ t ¶ x ¶ f (u) in which the characteristic velocity is and assumed to be positive, without loss ¶ u of generality. Here, the spatial discretization of Eq. (1) is given on an equally spaced one-dimensional mesh, in which the distance between two grid nodes are h, leading to an ordinary differential equation (ODE) system, i.e. du ¶ f = j = f ; i = 1; ; n: (2) x=x i dt ¶ x The numerical solutions of f will be discussed in the following sections, and the tem- poral solutions are given by using the third-order strongly stable Runge-Kutta method [28]. 2.1. The approximation of the ﬂux derivative The ﬁrst-order derivative of ﬂux function, i.e. f , can be approximated by using various numerical schemes. In a WCNS-type method, usually a central (compact) difference scheme is used to calculate the ﬂux derivatives at grid nodes, by utilizing the ﬂux function at midpoints and grid nodes. As aforementioned, because of their efﬁciency, in this work, explicit central dif- ference schemes which require nontridiagonal inversion are used to approximate the 4 derivative. Here, for achieving ﬁfth-order overall accuracy, three midpoint-to-node (MD) or midpoint-and-node-to-node (MND) formulas of sixth-order accuracy are given as a a a 1 2 3 e e e e e e f = ( f f )+ ( f f )+ ( f f ); (3) 1 1 3 3 5 5 i+ i i+ i i+ i 2 2 2 2 2 2 h h h b b b 1 2 3 e e e e e e f = ( f 1 f 1 )+ ( f f )+ ( f f ); (4) i+1 i1 i+2 i2 i+ i 2 2 h h h and c c c 1 2 3 e e e e e e f = ( f 1 f 1 )+ ( f f )+ ( f 3 f 3 ); (5) i+1 i1 i+ i i+ i h 2 2 h h 2 2 which are corresponding to the WCNS [29], Hybrid cell-edge and cell-node Weighted Compact Nonlinear Scheme (HWCNS) [30] and WCNS-midpoint-and-node-to-node difference (WCNS-MND) [16], respectively. Moreover, these three formulas show that they are not compact schemes, but can be taken as the special cases of the general WCNS [29]. The constant coefﬁcients in Eq. (3), (4) and (5), can be found in the corresponding references. It is worth noticed that the required midpoint stencils of HWCNS and WCNS- MND are more compact than that of Eq.(3), because they use node variables in the difference schemes, involving fewer midpoint interpolations. Moreover, Nonomura and Fujii [16] proved that WCNS-MND is more robust but less accurate than the orig- inal WCNS [5]and the original WCNS has higher propensity to blow up compared to WENO schemes when strong discontinuities are captured. WCNS-MND method is thus used in simulations of turbulence problems [31] and multi-species ﬂow problems [32]. In this work, most of the numerical results and discussions are given by using the sixth-order explicit MD scheme (Eq. (3)), and, of course, the other two formulas, Eq. (4) and (5), can be used straightforwardly. 2.2. Node-to-midpoint interpolation: WCNS As shown in the last subsection, midpoint ﬂux terms are unknown and should be evaluated before performing the MD or MND procedure. To achieve this goal, numer- 5 ical upwind ﬂux functions are used. Without loss of generality, a scalar form is used to explain the interpolation procedure. The scalar upwind ﬂux function is h i f = f (u )+ f (u ) ja ˆj u u ; (6) 1 1 1 1 1 i R;i L;i R;i L;i 2 2 2 2 2 where the subscript L and R indicate the variables at the left and right side of x , respectively, and a ˆ is the approximate eigenvalue. The question left to be settled is how to calculate u accurately, while maintaining numerical stability, nonoscillatory L=R;i and sharp discontinuity-capturing properties. For simplicity, we only consider the eval- uation of variables on the left side of x 1 , i.e. u 1 , in the following paragraphs. The i+ L;i+ 2 2 interpolations of the other midpoint variables are performed by using a symmetrical form of u 1 . L;i+ The ﬁfth-order node-to-midpoint reconstruction of WCNS to be used in this work was introduced in [5]. The basic idea of it is to approximate a linear optimal approxi- mation of the midpoint variable u 1 = u + (3u 20u 38u + 60u 5u ); (7) i i2 i1 i i+1 i+2 L;i+ 2 128 by using a nonlinear convex combination of lower-order interpolations. It is obvi- ous that a ﬁve-point full stencil S = fx ; x ; x ; x ; x g is used to calculate i2 i1 i i+1 i+2 i+ the high-order approximation. This optimal ﬁfth-order scheme can be equivalently represented by using three third-order polynomials each constructed on the following three-point substencil S =fx ; x ; x g; k = 1; 2; 3: (8) i+k3 i+k2 i+k1 i+ ;k Each of the third-order polynomials can be expressed in a generic form using the (approximated) nth derivatives (n = 1; 2) Dx (1) (2) u = u (x +Dx) = u + u Dx+ u ; (9) i i i i;k i;k L;i+ ;k where Dx = x x = . Speciﬁcally, the ﬁrst- and second-order derivatives are i+ 2 (1) u = (u 4u + 3u ); i2 i1 i i;1 2h (1) (10) u = (u u ); i+1 i1 i;2 2h (1) u = (3u + 4u u ); i i+1 i+2 i;3 2h 6 and (2) u = (u 2u + u ); i2 i1 i i;1 (2) (11) u = (u 2u + u ); i1 i i+1 i;2 (2) u = (u 2u + u ); i i+1 i+2 i;3 respectively. The linear optimal scheme is then represented as u = d u ; (12) 1 1 å k L;i+ L;i+ ;k 2 2 k=1 where the optimal linear weights are 1 10 5 d = ; d = ; d = : (13) 1 2 3 16 16 16 Obviously, directly using the optimal weights yields a ﬁfth-order accurate scheme. The resulting scheme, however, is inadequate to obtain a satisfying resolution for a steep gradient solution: it causes spurious oscillations. The nonlinear weights of Jiang and Shu [2] can be applied to replace the optimal linear weights in order to suppress non-physical oscillations. The so called JS nonlinear weights are deﬁned by a d k k w = ; a = ; (14) k k 3 2 (b +e) å k k=1 where the small parameter e = 10 is speciﬁed to avoid the denominator becoming zero, and b is the smoothness indicator, which is deﬁned as 2 2 (1) (2) b = hu + h u : (15) i;k i;k The main advantage of using nonlinear weights is that the contribution of oscillatory stencil(s) will be approximately eliminated in the ﬁnal interpolation while preserving those of relative smooth, and thus the spurious numerical oscillation can be suppressed. 2.3. Node-to-midpoint interpolation: ENO-like stencil-selection Instead of merely concentrating on improving the nonlinear weights, the ENO-like stencil-selection procedure [20] is introduced as an essential component of the pre- sented method. Firstly, the nonlinear smoothness measurement yielding strong scale- separation mechanism of the ﬁfth-order scheme is given as g = C + ; k = 1; 2; 3: (16) b +e 7 where t = jb b j is the global smooth indicator that was originally introduced in 5 1 3 [4], and the small threshold is deﬁned as e = 10 . Constant C = 1 is set, and the integer power q = 6 is used. As introduced by Fu et al. [20], larger integer power exponent q and smaller C are preferable for a stronger separation between resolved and non-resolved scales, and the discontinuity-detection capability can be signiﬁcantly enhanced. Rather than using the original nonlinear smoothness measurement, Eq. (16) is fur- ther normalized by c = ; (17) k=1 which is then subject to a sharp cutoff function 0; if c < C ; k T d = (18) 1; otherwise: The underlying idea of this design is that by introducing a value C as the threshold of smoothness, the candidate stencils are attributed as “smooth” or “oscillatory”, such that those stencils of oscillations are abandoned thereby, and only smooth ones are saved with their associated optimal linear weights used in the ﬁnal interpolation. The resulting weight functions are now given by d d (T) k k w = ; (19) d d k k k=1 where d denotes the optimal linear weight, and the above cut-off function d is incor- k k porated into the ﬁnal weighting to decide whether each candidate stencil is taken into account or not. Two promising advantages of this weighting procedure can be readily shown. Firstly, by removing the stencil crossing discontinuity completely, it ensures the numerical robustness of the scheme. Secondly, the background ﬁfth-order linear scheme can be fully recovered in smooth regions, including at smooth critical points. In fact, the stencil-selection procedure has simpliﬁed the nonlinear convex combi- nations of candidate stencils. Instead of using the continuous varying weights lead- ing to inﬁnite possible combinations, the presented method in fact uses only sev- eral candidate combinations, since each candidate stencil has only two possible re- 8 Table 1: The coefﬁcients of the equivalent single polynomial spatial reconstructions. if d = u ˆ S a a a a a 1;2;3 m;i2 m;i1 m;i m;i+1 m;i+2 1 m L;i+ ;m 1,1,1 u ˆ S 3/128 -5/32 45/64 15/32 -5/128 1 0 L;i+ ;0 0,1,1 u ˆ S 0 -1/12 5/8 1/2 -1/24 1 1 L;i+ ;1 1,1,0 u ˆ S 3/88 -5/22 75/88 15/44 0 1 2 L;i+ ;2 0,0,1 u ˆ S 0 0 3/8 3/4 -1/8 1 3 L;i+ ;3 0,1,0 u ˆ S 0 -1/8 3/4 3/8 0 1 4 L;i+ ;4 1,0,0 u ˆ S 3/8 -5/4 15/8 0 0 1 5 L;i+ ;5 1,0,1 u ˆ S 1/16 -5/24 5/8 5/8 -5/48 1 6 L;i+ ;6 Figure 1: Schematic of the equivalent candidate stencils of the ﬁfth-order nonlinear interpolation. The circle indicates that the stencil is not really continuous. sults of d d , i.e. d and 0. Therefore, similar to [24], each of the possible com- k k k binations of the nonlinear interpolation of u ˆ 1 can be represented as u ˆ = L;i+ L;i+ ;m a u + a u + a u + a u + a u , where the coefﬁcients are m;i2 i2 m;i1 i1 m;i i m;i+1 i+1 m;i+2 i+2 given in Table 1. The corresponding equivalent stencil of each convex combination is given in Fig.1, by which the meaning of ”stencil-selection” can be interpreted more clear. In addition, using the stencil-selection procedure results into several potential choices of high-order polynomials, and each can be optimized independently such as in [24]. Parameter C , serving as the global reference of smooth indicators, is an effective and a direct mean to control the spectral properties of the scheme for compressible turbu- lence simulation in which close embedded shocklets need to be captured [24]. Within 9 5 the scope of this work, however, the parameter is simply set as C = 10 for all the simulations without further investigation on the inﬂuence of its choice, but the advan- tages of importing the ENO-like stencil-selection procedure into compact nonlinear schemes are clearly showed. 3. Numerical results A variety of canonical problems are simulated to assess ﬁfth-order WCNS-JS, WCNS-Z, and the proposed scheme TCNS. One-dimensional linear advection equa- tion, one-dimensional inviscid Burgers equation and Euler equations of gas dynamics are used as model equations. The ideal-gas equation of state is given by p = (g 1)r e with g = 1:4 to close Euler equations. Node-to-midpoint interpolation is performed on characteristic variables to alleviate spurious oscillations [5]. Van leer scheme [8] is used for the computation of ﬂuxes. The proposed ﬁfth-order TCNS scheme is ﬁrst assessed through comparisons with the classical WCNS-JS, WCNS-Z on spectral prop- erties and ofﬁcial order of accuracy. Its rate of convergence is then veriﬁed by solv- ing one-dimensional linear advection equation, and one-dimensional inviscid Burgers equation. Performances of the above ﬁfth-order schemes in shock-capturing and wave resolutions are compared by solving selected test-problems both in one and two dimen- sions. The CFL= 0:6 is used as default for all numerical schemes and testing cases. 3.1. Approximate dispersion relation analysis Following the ADR analysis introduced by Pirozzoli [33] and Tu et al. [34], the spectral properties of different ﬁfth-order schemes are compared and shown in Fig. 2. The result of the presented method agrees very well with the background linear scheme for the low range of wave-numbers with a recovered wave-number reaching up to 1.76. A signiﬁcant improvement can also be found when cosmpared against the WCNS-JS in both dispersion and dissipation properties. While slight inferior result of dissipation property is given by TCNS when compared with WCNS-Z for the wave-number in intermediate range, TCNS performs better than WCNS-Z scheme for rest range of wave-numbers. Moreover, the performance of TCNS will be signiﬁcantly improved if 10 small C is used, and the adaptive dissipation control in [24] will further improve its overall performance. (a) (b) Figure 2: Approximated dispersion and dissipation properties of different ﬁfth-order schemes: dispersion (left) and dissipation (right). 3.2. Linear advection problem The one-dimensional Gaussian pulse advection problem [35] is used to assess the numerical order of accuracy of the proposed scheme. This problem is modeled by means of the linear advection equation, given by ¶ u ¶ u + = 0; x2 [0; 1]; (20) ¶ t ¶ x with periodic boundary-conditions and the initial conditions given by u(t; x = 0) = u(t; x = 1); (21) 300(xx ) u(0; x) = e ; x = 0:5: Time integration is performed up to t = 1, which corresponds to one period of the single wave propagation in time. A set of evenly distributed grids are progressively reﬁned by a factor of 2 from the most coarse grid with N = 51. Numerical simulation on each of the grid is conducted by using a time step size that even further reducing its value will not change the evaluated error of the numerical solution. 11 Table 2 and Fig. 3 illustrate the numerical error of those ﬁfth-order schemes used, where the result of TCNS scheme coincides with that of the corresponding linear scheme, indicating that it recovers the optimal order of convergence. A slight devi- ation can be found for WCNS-Z over the linear scheme for the coarse grids. WCNS-JS also shows approximate ﬁfth-order accuracy, but its resolution is signiﬁcantly lower than WCNS-Z and TCNS. In general, using the ENO-like stencil-selection procedure recovers the optimal linear scheme in this smooth ﬁeld. Table 2: L -error and convergence rate for different ﬁfth-order schemes used by solving the linear advection equation at t = 1 . Linear WCNS-JS WCNS-Z TCNS Error Order Error Order Error Order Error Order 51 5.22E-02 * 1.07E-01 * 5.67E-02 * 5.20E-02 * 101 3.30E-03 3.98 1.04E-02 3.37 3.57E-03 3.99 3.30E-03 3.98 201 1.16E-04 4.83 4.63E-04 4.49 1.20E-04 4.90 1.16E-04 4.83 401 3.69E-06 4.97 1.84E-05 4.66 3.72E-06 5.01 3.69E-06 4.97 801 1.16E-07 4.99 6.36E-07 4.85 1.16E-07 5.00 1.16E-07 4.99 1601 3.64E-09 4.99 2.02E-08 4.98 3.64E-09 4.99 3.64E-09 4.99 (a) L error (b) L error 2 ¥ Figure 3: Convergence rate of the L - and L -error for different ﬁfth-order schemes used by solving the 2 ¥ linear advection equation at t = 1. 12 3.3. Inviscid Burgers equation The one-dimensional inviscid Burgers equation [36, 37] is used to assess the actual order of accuracy of the proposed scheme when it is applied to a non-linear scalar equation. The governing equation is in the form of ¶ u ¶ 1 + u = 0; x2 [0; 2]; (22) ¶ t ¶ x 2 with periodic boundary-conditions and the initial conditions given by u(t; x = 0) = u(t; x = 1); (23) u(0; x) = + sin(p x): Exact solutions are computed by solving the derived general characteristic relation in reference [36]. The solutions are smooth for 0 t < 1=p , and a discontinuity devel- ops and starts to interact with the expansion wave if t 1=p . The results at t = 0:2 and t = 0:7 are both given to show the continuous and discontinuous distributions. Solu- tions obtained by the above ﬁfth-order schemes are compared against the exact solution in Fig. 4, in which, good agreement is presented both for the smooth and the discon- tinuous solution. Moreover, at t = 0:2, L -error and convergence rate for each scheme are presented in Table 3. All of the schemes can achieve ﬁfth-order of accuracy as the grids are reﬁned. Perfect agreement can be found between the TCNS and the back- ground ﬁfth-order linear scheme. WCNS-Z scheme shows minor deviation compared with TCNS scheme only when the grid is relatively coarse. WCNS-JS scheme is less accurate than the other ﬁfth-order schemes. Again, TCNS scheme recovers the optimal linear scheme if the solution is smooth. Moreover, it should be noticed that in the discontinuous solution, minor overshoot has been found. This problem can be avoided by using the MND method, as shown in Fig.4(b). At the meantime, MND method reduces the resolution of the result, as already introduced by Nonomura and Fujii [16]. Whereas, in this work, only the node- to-midpoint interpolation methods are discussed, and thus the detail of MND method is not further investigated. As mentioned above, the presented method can be applied with using all the midpoint-to-node or midpoint-and-node-to-node methods. 13 (a) t = 0:2 (b) t = 0:7 Figure 4: One-dimensional inviscid Burgers equation: numerical solutions and exact solution at t = 0:2 (left) and t = 0:7 (right). Grid resolution 161. Table 3: L -error and convergence rate for different ﬁfth-order schemes used by solving the 1-D inviscid Burgers equation at t = 0.2. Linear WCNS-JS WCNS-Z TCNS Error Order Error Order Error Order Error Order 41 9.96E-04 * 1.41E-03 * 1.01E-03 * 9.96E-04 * 81 9.04E-05 3.46 1.34E-04 3.40 9.07E-05 3.48 9.04E-05 3.46 161 3.21E-06 4.81 5.06E-06 4.73 3.22E-06 4.82 3.21E-06 4.81 321 1.00E-07 5.00 1.63E-07 4.96 1.00E-07 5.00 1.00E-07 5.00 641 2.82E-09 5.15 4.82E-09 5.07 2.82E-09 5.15 2.82E-09 5.15 3.4. Sod and Lax shock tube problem Riemann initial value problems of Sod [38] and Lax [39] are used to evaluate the shock-capturing capability of the proposed scheme. The initial conditions for the Sod problem is (1; 0; 1) x2 [0; 0:5]; (r; u; p) = (24) (0:125; 0; 0:1) x2 (0:5; 1]: and the ﬁnal time is t = 2 by solving the problem on an evenly distributed grid of N = 101 points. 14 (a) (b) Figure 5: Sod problem: numerical solutions and the exact solution at t = 0:2. Overall, good agreement of numerical solutions with the exact solution is achieved for all of the schemes shown in Fig. 5(a). Steep gradient solutions such as shock wave and contact discontinuity are captured without oscillations. In Fig. 5(b), WCNS-Z and TCNS schemes yield a relative sharp proﬁle compared to WCNS-JS in the vicinity of contact discontinuity, and TCNS is found with a better performance of resolving the contact discontinuity when compared with WCNS-Z. The initial conditions for the Lax shock-tube problem is (0:445; 0:698; 3:528) x2 [0; 0:5]; (r; u; p) = (25) (0:5; 0; 0:571) x2 (0:5; 1]: This case is run on a grid of N = 101 points with uniform distribution, and the results at t = 0:14 are given. Density and velocity distributions are compared in Fig.6(a), and Fig. 6(b) respec- tively. WCNS-Z overshoots the velocity slightly at the tail of the expansion fan, but TCNS shows accurate and oscillation-free result. WCNS-JS smears the shock wave signiﬁcantly, indicating extra numerical dissipation. 15 (a) Density (b) Velocity Figure 6: Lax problem: numerical solutions and the exact solution at t = 0:14. 3.5. Shock-density wave interaction The shock-density wave interaction problem [40] is characterized by a right mov- ing Mach 3 shock interacting with sine waves in density ﬁeld. The multi-scale wave structure is evolved after the shock wave interacts with the oscillating density wave, and both the shock-capturing and wave-resolution capabilities are evaluated thereafter. The problem is initialized by (3:857; 2:629; 10:333); x2 [0; 1]; (r; u; p) = (26) (1+ 0:2sin(5x); 0; 1); x2 (1; 10]: This case is run on a grid of N = 201 points which are uniformly distributed and the ﬁnal time is t = 1:8. Numerical solution of WCNS-JS on a grid of N = 2001 is used as the reference ”exact” solution. As shown in Fig. 7, TCNS produces considerably better resolved density waves be- hind the shock wave compared with WCNS-JS and WCNS-Z. Particularly, the result of WCNS-JS indicates strong numerical dissipation, since the small-scale wave structure is signiﬁcantly smeared. Moreover, TCNS does not produce spurious oscillations in the ﬂow ﬁeld. 16 (a) (b) Figure 7: Shock-density wave interaction problem: numerical solutions and the exact solution at t = 1:8. 3.6. Two-dimensional Riemann problem In this section, the conﬁgurations 3 and 6 out of 19 2-D Riemann problems used by Lax and Liu [41] are speciﬁcally taken to evaluate performances of the proposed numerical schemes. 3.6.1. The conﬁguration 3 The initial conditions for the conﬁguration 3 are given by 1 1 (1:5; 0; 0; 1:5) x2 ; 1 and y2 ; 1 ; > 2 2 1 1 (0:5323; 1:206; 0; 0:3) x2 0; and y2 ; 1 ; 2 2 (r; u; v; p) = (27) > 1 1 >(0:138; 1:206; 1:206; 0:029) x2 0; and y2 0; ; 2 2 : 1 1 (0:5323; 0; 1:206; 0:3) x2 ; 1 and y2 0; : 2 2 Boundary conditions are given by ¶u(t; x; y) = 0; x = 0; 1; 8t; y; ¶ x (28) ¶u(t; x; y) = 0; y = 0; 1; 8t; x: ¶ y Figure 8 presents the results obtained by using WCNS-JS, WCNS-Z and TCNS schemes, at t = 0:3. The density proﬁles are best predicted by TCNS. WCNS-Z yields 17 a solution very close to that of TCNS, but TCNS captures more vortexes along the contact lines. WCNS-JS smears a few small scales that can be resolved by TCNS and WCNS-Z due to relative large dissipations. (a) WCNS-JS (b) WCNS-Z (c) TCNS Figure 8: Conﬁgurations 3 of 2-D Riemann problems in [41]: 30 density contours ranging from 0.1 to 1.8 at t = 0:23 obtained on a grid of 1024 2014. 18 3.6.2. The conﬁguration 6 The initial conditions for the conﬁguration 6 are given by 1 1 (1; 0:75;0:5; 1) x2 ; 1 and y2 ; 1 ; > 2 2 1 1 (2; 0:75; 0:5; 1) x2 0; and y2 ; 1 ; 2 2 (r; u; v; p) = (29) 1 1 > (1;0:75; 0:5; 1) x2 0; and y2 0; ; 2 2 : 1 1 (3;0:75;0:5; 1) x2 ; 1 and y2 0; : 2 2 Boundary conditions are the same as in the last case. Results obtained by using WCNS-JS, WCNS-Z and TCNS schemes at t = 0:3 are shown in Fig. 9. The results of WCNS-JS and WCNS-Z are similar. Whereas, using TCNS obtains abundant small scale ﬂow structures along the contact lines, indicating lower numerical dissipation. 19 (a) WCNS-JS (b) WCNS-Z (c) TCNS Figure 9: Conﬁgurations 6 of 2-D Riemann problems in [41]: 40 density contours ranging from 0.1 to 2.9 at t = 0:3 obtained on a grid of 1024 2014. 3.7. Shock vortex interaction This problem was modeled in two dimensions by Jiang and Shu [2]. It involves a vortex perturbing a stationary shock. The computational domain is [0; 2] [0; 1]. This ﬁeld is discretized by 251 101 grid. Initially, a stationary shock is located at x = 0:5 normal to the x direction. The left state of this shock is speciﬁed as (r; u; v; p) = (1; 1:1 g; 0; 1). A small vortex centered at (0.25, 0.5) is superposed to the ﬂowﬁled on the left hand side of the normal shock. The superposition is performed through 20 perturbations on velocity (u; v), temperature T , given by T = p=r , and entropy S, deﬁned as S = ln( p=r ) of the mean ﬂow. Speciﬁcally, the perturbation variables are a(1t ) u ˜ = et e sinq a(1t ) u ˜ =et e cosq (30) 2 2a(1g ) (g1)e e > T = 4ag S = 0 2 2 where t = r=r , and r = (x x ) +(y y ) , r = 0:05, e = 0:3, a = 0:204 are taken c c c c from the reference [2]. The results at t = 0:6 are presented in Fig.10. (a) WCNS-JS (b) WCNS-Z (c) TCNS Figure 10: Shock vortex interaction problem: 90 density contours ranging from 1.19 to 1.37. 21 The ﬂow contours of different schemes are similar. The density distribution along the central line of the ﬂow ﬁeld are further shown in Fig.11, where the result obtained by using WCNS-Z on a reﬁned mesh of 1001 401 is used as the “exact” solution. It can be found that the captured vortex of TCNS is more accurate, comparing with the other two results of WCNS-JS and WCNS-Z. Figure 11: Shock vortex interaction problem: density distribution along y = 0:5. 3.8. Rayleigh-Taylor instability Rayleigh-Taylor instability problem which contains both discontinuities and com- plex ﬂow structures [42], is also used to examine the numerical dissipation of the pre- sented method. The initial conditions are given by (2; 0;0:025a cos(8p x); 1+ 2y) x2 [0; 0:25] and y2 [0; 0:5); (r; u; v; p) = (1; 0;0:025a cos(8p x); 1+ 3=2) x2 [0; 0:25] and y2 [0:5; 1]; (31) where a is the speed of sound, given by a = g and a different g = is used for r 3 this speciﬁc case. Reﬂecting boundary conditions are imposed at the left and right side of the domain, and constant boundary conditions are given for the top and the bottom 22 sides, in details (1; 0; 0; 2:5) y = 1; 8t; x; (r; u; v; p) = (32) (2; 0; 0; 1) y = 0; 8t; x: Two source terms r , and r v are added to the right hand side of the third and the fourth equation, respectively. Two sets of grids are used, i.e. 128 512 and 256 1024. Density proﬁles at ﬁnal time t = 1:95 calculated by WCNS-JS, WCNS-Z and TCNS are shown in Fig. 12. 23 (a) Grid resolution 128 512 (b) Grid resolution 256 1024 Figure 12: Rayleigh-Taylor instability problem: 30 density contours ranging from 0.9 to 2.2. 24 TCNS has resolved much more abundant vortical structures than WCNS-JS and WCNS-Z, suggesting that TCNS is signiﬁcantly less dissipative. Furthermore, the low-dissipation property of TCNS induces symmetry breaking of the ﬂow ﬁeld. In par- ticular, the small structures resolved by using TCNS with grid resolution of 128 512 are comparable to that of using WCNS-JS with grid resolution of 256 1024. Even the result of WCNS-Z only shows little advantage (if any) while doubling the grid resolution. 4. Conclusions As a subsequent work of compact nonlinear scheme and weighed compact non- linear scheme, we have introduced a novel compact nonlinear scheme in this article, by applying the ENO-like stencil-selection procedure. The method is targeted optimal linear interpolation in the node-to-midpoint interpolation step of the compact nonlin- ear scheme, and thus the method is named as targeted CNS (TCNS). Numerical results of one-dimensional scalar equations and one- or two-dimensional Euler equations are given to examine the performance of the method, involving strong discontinuities and broadband ﬂuctuations. The presented method uses a novel smoothness measurement which possesses strong scale separation mechanism, and the ENO-like stencil-selection procedure is capable to target optimal linear weights while maintaining excellent shock-capturing capability. ADR analysis shows that TCNS recovers the background linear scheme up to higher wave number, even using a relative larger C . Examining the numerical results, it can be found that TCNS has recovered optimal linear scheme in smooth ﬁeld, where the optimal linear weights d is directly applied in the node-to-midpoint interpolation. Whereas, WCNSs, including WCNS-JS and WCNS-Z, approximate the optimal linear weights but can not perfectly recover them. Therefore, the presented TCNS is capa- ble to capture much more abundant wave structures in the simulations. Especially in the Rayleigh-Taylor instability problem, TCNS achieves similar or even better result comparing with WCNS-JS on a coarser grid. 25 Acknowledgments Financial support for this work was provided through grant number 20187413071020008. References [1] X. D. Liu, S. Osher, T. Chan, Weighted essentially non-oscillatory schemes, Jour- nal of Computational Physics 115 (1) (1994) 200–212. [2] G.-S. Jiang, C.-W. Shu, Efﬁcient implementation of weighted ENO schemes, Journal of Computational Physics 126 (1) (1996) 202–228. [3] A. K. Henrick, T. D. Aslam, J. M. Powers, Mapped weighted essentially non- oscillatory schemes: Achieving optimal order near critical points, Journal of Computational Physics 207 (2) (2005) 542–567. [4] R. Borges, M. Carmona, B. Costa, W. S. Don, An improved weighted essentially non-oscillatory scheme for hyperbolic conservation laws, Journal of Computa- tional Physics 227 (6) (2008) 3191–3211. doi:10.1016/j.jcp.2007.11.038. [5] X. Deng, H. Zhang, Developing high-order weighted compact nonlinear schemes, Journal of Computational Physics 165 (1) (2000) 22 – 44. doi:10.1006/jcph. 2000.6594. [6] S. K. Lele, Compact ﬁnite difference schemes with spectral-like resolution, Journal of Computational Physics 103 (1) (1992) 16 – 42. doi:10.1016/ 0021-9991(92)90324-R. [7] P. L. Roe, Approximate Riemann solvers, parameter vectors, and difference schemes, Journal of Computational Physics 43 (2) (1981) 357–372. [8] B. van Leer, Flux-vector splitting for the Euler equations, in: E. Krause (Ed.), Eighth International Conference on Numerical Methods in Fluid Dynamics, Vol. 170 of Lecture Notes in Physics, Springer Berlin Heidelberg, Aachen, Germany, 1982, pp. 507–512. 26 [9] M. S. Liou, On a new class of ﬂux splittings, in: M. Napolitano, F. Sabetta (Eds.), Thirteenth International Conference on Numerical Methods in Fluid Dynamics, Vol. 414 of Lecture Notes in Physics, Springer Berlin Heidelberg, Berlin, Heidel- berg, 1993, pp. 115–119. [10] T. Nonomura, N. Iizuka, K. Fujii, Freestream and vortex preservation properties of high-order weno and wcns on curvilinear grids, Computers & Fluids 39 (2) (2010) 197 – 214. doi:10.1016/j.compfluid.2009.08.005. [11] D. Wang, X. Deng, G. Wang, Y. Dong, Developing a hybrid ﬂux function suitable for hypersonic ﬂow simulation with high-order methods, Int. J. Numer. Meth. Fluids 81 (5) (2016) 309–327. [12] G.-H. TU, J.-Q. CHEN, M.-L. MAO, X.-H. ZHAO, H.-Y. LIU, On the splitting methods of inviscid ﬂuxes for implementing high-order weighted compact non- linear schemes, Applied Mathematics and Mechanics 37 (12) (2016) 1324–1344. doi:10.21656/1000-0887.370518. [13] H. Zhu, Z. Yan, H. Liu, M. Mao, Properties of osher ﬂux with entropy ﬁx in high- order WCNS, Acta Aeronautica et Astronautica Sinica 38 (5) (2017) 120543. [14] X. Deng, X. Liu, M. Mao, H. Zhang, Investigation on weighted compact ﬁfth- order nonlinear scheme and applications to complex ﬂow, in: 17th AIAA Com- putational Fluid Dynamics Conference, no. 0 in Fluid Dynamics and Co-located Conferences, AIAA Paper 2005-5246, 2005. doi:10.2514/6.2005-5246. [15] T. Nonomura, K. Fujii, Effects of difference scheme type in high-order weighted compact nonlinear schemes, Journal of Computational Physics 228 (10) (2009) 3533 – 3539. doi:10.1016/j.jcp.2009.02.018. [16] T. Nonomura, K. Fujii, Robust explicit formulation of weighted compact non- linear scheme, Computers & Fluids 85 (2013) 8 – 18. doi:10.1016/j. compfluid.2012.09.001. 27 [17] Z. Yan, H. Liu, M. Mao, H. Zhu, X. Deng, New nonlinear weights for improving accuracy and resolution of weighted compact nonlinear scheme, Computers & Fluids 127 (2016) 226 – 240. doi:10.1016/j.compfluid.2016.01.005. [18] Z.-G. Yan, H. Liu, Y. Ma, M. Mao, X. Deng, Further improvement of weighted compact nonlinear scheme using compact nonlinear interpolation, Computers & Fluids 156 (2017) 135 – 145. doi:10.1016/j.compfluid.2017.06.028. [19] Y. Ma, H. Liu, Z. Yan, M. Mao, X. Deng, Investigation of compact interpolation based on hwcns schemes, Chinese Journal of Computational Mechanics 32 (3) (2015) 388–393. [20] L. Fu, X. Y. Hu, N. A. Adams, A family of high-order targeted ENO schemes for compressible-ﬂuid simulations, Journal of Computational Physics 305 (2016) 333–359. [21] L. Fu, X. Y. Hu, N. A. Adams, A targeted eno scheme as implicit model for tur- bulent and genuine subgrid scales, Journal of Computational Physics 349 (2017) 97–121. [22] O. Haimovich, S. H. Frankel, Numerical simulations of compressible multicom- ponent and multiphase ﬂow using a high-order targeted ENO (TENO) ﬁnite- volume method, Computers & Fluids 146 (2017) 105–116. doi:10.1016/j. compfluid.2017.01.012. [23] L. Fu, X. Hu, N. A. Adams, A targeted eno scheme as implicit model for turbulent and genuine subgrid scales, Communications in Computational Physics. [24] L. Fu, X. Y. Hu, N. A. Adams, A new class of adaptive high-order tar- geted ENO schemes for hyperbolic conservation laws, Journal of Computational Physicsdoi:10.1016/j.jcp.2018.07.043. [25] Z. Sun, S. Inaba, F. Xiao, Boundary variation diminishing (bvd) reconstruction: A new approach to improve godunov schemes, Journal of Computational Physics 322 (2016) 309 – 325. doi:10.1016/j.jcp.2016.06.051. 28 [26] X. Deng, H. Maekawa, Compact high-order accurate nonlinear schemes, Journal of Computational Physics 130 (1) (1997) 77 – 91. doi:10.1006/jcph.1996. [27] S. Pirozzoli, On the spectral properties of shock-capturing schemes, Journal of Computational Physics 219 (2) (2006) 489 – 497. doi:10.1016/j.jcp.2006. 07.009. [28] S. Gottlieb, C. Shu, E. Tadmor, Strong stability-preserving high- order time discretization methods, SIAM Review 43 (1) (2001) 89– 112. arXiv:https://doi.org/10.1137/S003614450036757X, doi:10.1137/S003614450036757X. URL https://doi.org/10.1137/S003614450036757X [29] X. Deng, High-order accurate dissipative weighted compact nonlinear schemes, Science in China Series A: Mathematics 45 (3) (2002) 356. doi:10.1360/ 02ys9037. [30] X. Deng, New high-order hybrid cell-edge and cell-node weighted compact non- linear schemes, in: 20th AIAA Computational Fluid Dynamics Conference, AIAA Paper 2011-3857, 2011. doi:10.2514/6.2011-3857. [31] G. Zhao, M. Sun, S. Xie, H. Wang, Numerical dissipation control in an adaptive wcns with a new smoothness indicator, Applied Mathematics and Computation 330 (2018) 239 – 253. doi:10.1016/j.amc.2018.01.019. [32] M. L. Wong, S. K. Lele, High-order localized dissipation weighted compact non- linear scheme for shock- and interface-capturing in compressible ﬂows, Journal of Computational Physics 339 (2017) 179 – 209. doi:10.1016/j.jcp.2017. 03.008. [33] S. Pirozzoli, On the spectral properties of shock-capturing schemes, Journal of Computational Physics 219 (2006) 489–497. doi:10.1016/j.jcp.2006.07. 29 [34] G.-H. Tu, X.-G. Deng, M.-L. Mao, Spectral property comparison of ﬁfth-order nonlinear wcns and weno difference schemes, Acta Aerodynamica Sinica 30 (6) (2013) 709–712. [35] N. K. Yamaleev, M. H. Carpenter, A systematic methodology for constructing high-order energy stable weno schemes, Journal of Computational Physics 228 (11) (2009) 4248 – 4272. doi:https://doi.org/10.1016/j.jcp. 2009.03.002. URL http://www.sciencedirect.com/science/article/pii/ S0021999109001132 [36] A. Harten, B. Engquist, S. Osher, S. R. Chakravarthy, Uniformly high order ac- curate essentially non-oscillatory schemes, iii, Journal of Computational Physics 131 (1) (1997) 3 – 47. doi:https://doi.org/10.1006/jcph.1996.5632. URL http://www.sciencedirect.com/science/article/pii/ S0021999196956326 [37] G. Gerolymos, D. Snchal, I. Vallet, Very-high-order weno schemes, Journal of Computational Physics 228 (23) (2009) 8481 – 8524. doi:https://doi.org/10.1016/j.jcp.2009.07.039. URL http://www.sciencedirect.com/science/article/pii/ S0021999109003908 [38] G. A. Sod, A survey of several ﬁnite difference methods for systems of nonlinear hyperbolic conservation laws, Journal of Computational Physics 27 (1) (1978) 1 – 31. doi:https://doi.org/10.1016/0021-9991(78)90023-2. URL http://www.sciencedirect.com/science/article/pii/ [39] P. D. Lax, Weak solutions of nonlinear hyperbolic equations and their numerical computation, Communications on Pure and Applied Mathematics 7 (1) 159–193. arXiv:https://onlinelibrary.wiley.com/doi/pdf/10.1002/cpa. 3160070112, doi:10.1002/cpa.3160070112. 30 URL https://onlinelibrary.wiley.com/doi/abs/10.1002/cpa. [40] C.-W. Shu, S. Osher, Efﬁcient implementation of essentially non-oscillatory shock-capturing schemes, ii, Journal of Computational Physics 83 (1) (1989) 32 – 78. doi:https://doi.org/10.1016/0021-9991(89)90222-2. URL http://www.sciencedirect.com/science/article/pii/ [41] P. Lax, X. Liu, Solution of two-dimensional riemann problems of gas dynam- ics by positive schemes, SIAM Journal on Scientiﬁc Computing 19 (2) (1998) 319–340. arXiv:https://doi.org/10.1137/S1064827595291819, doi: 10.1137/S1064827595291819. URL https://doi.org/10.1137/S1064827595291819 [42] Z. Xu, C.-W. Shu, Anti-diffusive ﬂux corrections for high order ﬁnite difference weno schemes, Journal of Computational Physics 205 (2) (2005) 458 – 485. doi:https://doi.org/10.1016/j.jcp.2004.11.014. URL http://www.sciencedirect.com/science/article/pii/ S0021999104004760

Journal

Physics – arXiv (Cornell University)

Published: Oct 16, 2018

References

Weighted essentially non-oscillatory schemes

Liu
Efficient implementation of weighted ENO schemes

Jiang
Mapped weighted essentially non-oscillatory schemes: achieving optimal order near critical points

Henrick
An improved weighted essentially non-oscillatory scheme for hyperbolic conservation laws

Borges
Developing high-order weighted compact nonlinear schemes

Deng
A high-resolution hybrid compact-ENO scheme for shock-turbulence interaction problems

Adams
An efficient low-dissipation hybrid weighted essentially non-oscillatory scheme

Hu
An efficient adaptive high-order scheme based on the WENO process

Sun
Compact finite difference schemes with spectral-like resolution

Lele
Approximate Riemann solvers, parameter vectors, and difference schemes

Roe
Flux-vector splitting for the Euler equations

van Leer
On a new class of flux splittings

Liou
Freestream and vortex preservation properties of high-order WENO and WCNS on curvilinear grids

Nonomura
Developing a hybrid flux function suitable for hypersonic flow simulation with high-order methods

Wang
On the splitting methods of inviscid fluxes for implementing high-order weighted compact nonlinear schemes

Tu
Properties of Osher flux with entropy fix in high-order WCNS

Zhu
Investigation on weighted compact fifth-order nonlinear scheme and applications to complex flow

Deng
Effects of difference scheme type in high-order weighted compact nonlinear schemes

Nonomura
Robust explicit formulation of weighted compact nonlinear scheme

Nonomura
New nonlinear weights for improving accuracy and resolution of weighted compact nonlinear scheme

Yan
Further improvement of weighted compact nonlinear scheme using compact nonlinear interpolation

Yan
Investigation of compact interpolation based on HWCNS schemes

Ma
A family of high-order targeted ENO schemes for compressible-fluid simulations

Fu
Numerical simulations of compressible multicomponent and multiphase flow using a high-order targeted ENO (TENO) finite-volume method

Haimovich
A simple extending strategy for TENO scheme: scalar equation and Euler equations

Zhang
Detonation simulations with a fifth-order TENO scheme

Dong
Boundary variation diminishing (BVD) reconstruction: a new approach to improve Godunov schemes

Sun
On the spectral properties of shock-capturing schemes

Pirozzoli
Strong stability-preserving high-order time discretization methods

Gottlieb
High-order accurate dissipative weighted compact nonlinear schemes

Deng
New high-order hybrid cell-edge and cell-node weighted compact nonlinear schemes

Deng
Numerical dissipation control in an adaptive WCNS with a new smoothness indicator

Zhao
High-order localized dissipation weighted compact nonlinear scheme for shock- and interface-capturing in compressible flows

Wong
High order weighted essentially non-oscillatory WENO-Z schemes for hyperbolic conservation laws

Castro
A new class of adaptive high-order targeted ENO schemes for hyperbolic conservation laws

Fu
Study of quasi-linear spectral analysis method of high-order weighted nonlinear schemes

Mao
A systematic methodology for constructing high-order energy stable WENO schemes

Yamaleev
Uniformly high order accurate essentially non-oscillatory schemes, III

Harten
Very-high-order WENO schemes

Gerolymos
A survey of several finite difference methods for systems of nonlinear hyperbolic conservation laws

Sod
Weak solutions of nonlinear hyperbolic equations and their numerical computation

Lax
Efficient implementation of essentially non-oscillatory shock-capturing schemes, II

Shu
Solution of two-dimensional Riemann problems of gas dynamics by positive schemes

Lax
Anti-diffusive flux corrections for high order finite difference WENO schemes

Xu

Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Towards optimal high-order compact schemes for simulating compressible flows

Towards optimal high-order compact schemes for simulating compressible flows

Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Towards optimal high-order compact schemes for simulating compressible flows

Towards optimal high-order compact schemes for simulating compressible flows

References

Abstract

Journal

Recommended Articles

References

Our policy towards the use of cookies