Two-Stage Fourth-order Gas-kinetic Scheme for Three-dimensional Euler and Navier-Stokes Solutions
For the one-stage third-order gas-kinetic scheme (GKS), success applications have been achieved for the three-dimensional compressible flow computations GKS-high-1 (). The high-order accuracy of the scheme is obtained directly by integrating a multidimensional time-accurate gas distribution function over the cell interface within a time step without implementing Gaussian quadrature points and Runge-Kutta time-stepping technique. However, for the further increasing the order of the scheme, such as the fourth-order one, the formulation becomes very complicated for the multidimensional flow. Recently, a two-stage fourth-order GKS with high efficiency has been constructed for two-dimensional inviscid and viscous flow computations GRP-high (); GKS-high (), and the scheme uses the time accurate flux function and its time derivatives. In this paper, a fourth-order GKS is developed for the three-dimensional flows under the two-stage framework. Based on the three-dimensional WENO reconstruction and flux evaluation at Gaussian quadrature points on a cell interface, the high-order accuracy in space is achieved first. Then, the two-stage time stepping method provides the high accuracy in time. In comparison with the formal third-order GKS GKS-high-1 (), the current fourth-order method not only improves the accuracy of the scheme, but also reduces the complexity of the gas-kinetic solver greatly. More importantly, the fourth-order GKS has the same robustness as the second-order shock capturing scheme. This scheme is applied to both inviscid and viscous, and low and high speed flow computations. Numerical results validate the outstanding reliability and applicability of the scheme for three-dimensional flows, such as turbulent one. With the count of the degree of freedom in simulating a flow field, the current two-stage fourth-order multidimensional GKS becomes one of the most reliable and efficient higher-order schemes among all shock capturing higher-order schemes for computational fluid dynamics.
keywords:two-stage fourth-order discretization, gas-kinetic scheme, Navier-Stokes solutions, shock capturing scheme
Over the past half century, the computational fluid dynamics has been developed into a powerful tool for solving the fluid flow problems in industrial applications. Currently, there are a gigantic amount of numerical methods in literature. In comparison with well-developed second-order schemes, the higher-order methods can provide more accurate solutions, but they are less robust and more complicated. For the high-order schemes, there are generally three parts, i.e. spatial discretization, temporal discretization, and flux solvers. For the spatial discretization, many high-order methods have been developed, including the discontinuous Galerkin (DG) DG-1 (); DG-2 (); DG-3 (), essential non-oscillatory (ENO) ENO-1 (); ENO-2 (), weighted essential non-oscillatory (WENO) WENO (); WENO-JS (); WENO-Z (), and many others. For most of those methods, the exact and approximate Riemann solvers Riemann-appro () are used for flux evaluation. Due to the first-order evolution of Riemann solvers, the Runge-Kutta method is implemented to achieve higher order temporal accuracy TVD-rk (), in which -stage is needed for -th order accuracy. Instead of Riemann solver, many schemes have been developed based on the time-dependent flux function with high-order evolution, such as the generalized Riemann problem solvers GRP1 (); GRP2 (); GRP3 (); GRP-high-2 (), gas-kinetic scheme GKS-Xu1 (); GKS-Xu2 (), and AEDR methods Riemann-appro (); ADER (). High-order temporal accuracy can be achieved in a one-stage framework with the time integration of the time-dependent flux function. Recently, in order to increase the accuracy and efficiency of these schemes a two-stage fourth-order method has been developed for the time-dependent flux solvers GRP-high (); GKS-high (), where both the flux and temporal derivative of flux function are used in the construction of higher-order schemes. By combining the multi-stage multi-derivative technique, a family of higher-order schemes has been constructed as well MSMD-GKS ().
In the past decades, the gas-kinetic scheme (GKS) based on the Bhatnagar-Gross-Krook (BGK) model BGK-1 (); BGK-2 (); BGK-3 () has been developed systematically for the compressible flow computations GKS-Xu1 (); GKS-Xu2 (); GKS-Xu3 (). Different from the traditional CFD methods based on the macroscopic governing equations, the main advantages of the gas-kinetic scheme are the followings. (i) The inviscid and viscous coupling in the flux evolution GKS-Xu1 (); GKS-Xu2 (); (ii) Multi-dimensionality with the inclusion of both normal and tangential gradients of flow variables in the flux function across a cell interface GKS-high-4 (); (iii) Compact stencils can be constructed with the use of the time accurate cell interface flow variables at the next time level GKS-high-4 (); (iv) Extension to the whole flow regime from the rarefied to the continuum one UGKS-1 (); UGKS-2 (). Recently, with the high-order initial reconstruction, the third-order gas-kinetic schemes have been constructed GKS-high-2 (); GKS-high-3 (); GKS-high-4 (), in which the flux evaluation is based on the moments of spatial-temporal dependent gas distribution function. High-order accuracy can be achieved in a one-stage scheme without Gaussian point integration for spatial accuracy and Runge-Kutta technique for temporal accuracy. However, with the one-stage gas evolution model, the formulation of gas-kinetic scheme can become very complicated for the further development, such as the fourth-order method GKS-high-5 (), especially for multidimensional computations. The two-stage fourth-order temporal discretization for time-accurate flux solvers in GRP-high () provides a reliable framework to further develop the GKS into fourth-order and even higher accuracy with the implementation of the traditional second-order or third-order flux functions GKS-high (); GKS-high-6 (); MSMD-GKS (). Most importantly, this scheme is robust, and works perfectly from the subsonic to the hypersonic viscous heat conducting flows. The robustness of the GKS in comparison with Riemann solver based CFD methods is solely due to the differences in the dynamical evolution model of the flux function. The GKS flux follows the dynamics from the particle free transport, to including collisional effect, and to the NS distribution function with the count of intensive particle collisions as the ratio of the time step over the particle collision time becoming large. The real physics used in the flux depends on the local xu-liu (). In real NS computations for the compressible high Mach number flow, the ratio of is not too large as people think of for the validation of the NS modeling liu-kn (). However, for the Riemann solver based CFD methods, at the beginning of the step it is already assumed that there are infinity number of particle collision to generate distinguishable waves in the Riemann solution, and the collision needs to be reduced for the NS solutions. Theoretically, there is no such a physical process for the Riemann solver-based Godunov-type schemes for the NS equations. For the second-order schemes, it is hard to distinguish the dynamical differences from the GKS and Riemann solver. However, for the higher-order schemes it seems that a reliable physical evolution model becomes more important due to the absence of large numerical dissipation in the second-order schemes, and the delicate flow structures captured in higher-order schemes depend on the quality of the solvers greatly MSMD-GKS ().
In this paper, with the two-stage fourth-order discretization, a multidimensional fourth-order gas-kinetic scheme is constructed for simulating three-dimensional flows. High-order accuracy in space is achieved by the three-dimensional WENO method WENO (); WENO-JS (); WENO-Z () and Gaussian quadrature points for the numerical fluxes. In comparison with the formal three-dimensional scheme GKS-high-1 (), the current fourth-order scheme reduces the complexity of the gas-kinetic flux solver greatly, and improves the robustness of scheme. Many numerical tests, including both inviscid and viscous, and low and high speed flow computations, will be used to validate the current fourth-order method. Numerical results show that the current scheme has the same reliability and applicability as the well-developed second-order scheme, but is much more accurate and effective to capture the complicated flow structures. With the count of the degrees of freedom for the description of a flow field, the current scheme provides state-of-art solutions from a higher-order scheme from the incompressible to the hypersonic flow simulations.
This paper is organized as follows. In Section 2, BGK equation and finite volume scheme are briefly reviewed. The general formulation for the two-stage temporal discretization is introduced in Section 3, and the procedure of spatial reconstruction is given in Section 4. Section 5 includes numerical examples to validate the current algorithm. The last section is the conclusion.
2 BGK equation and finite volume scheme
where is the gas distribution function, is the corresponding equilibrium state, and is the collision time. The collision term satisfies the compatibility condition
where , the internal variables is equal to , , is the degrees of freedom, i.e. for three-dimensional flows and is the specific heat ratio. In the continuum region, the gas distribution function can be expanded as
where . Based on the Chapman-Enskog expansion, the macroscopic equations can be derived GKS-Xu1 (); GKS-Xu2 (). For the Euler equations, the zeroth-order truncation is taken, i.e. . For the Navier-Stokes equations, the first-order truncation is
With the higher order truncations, the Burnett and super-Burnett equations can be derived.
Taking moments of Eq.(1) and integrating over the control volume with , the semi-discretized form of finite volume scheme can be written as
where are the conservative flow variables, is the cell averaged value over control volume and . For the three-dimensional computation, the Gaussian quadrature for the numerical fluxes is used to achieve the accuracy in space, and the numerical fluxes in -direction is given as an example
where , are the Gauss quadrature points for and are corresponding weights, and is used in this paper. Based on the spatial reconstruction, which will be presented in the following section, the reconstructed point value and the spatial derivatives at each Gauss quadrature points can be obtained and the numerical fluxes can be provided by the flow solvers. Similarly, the numerical fluxes in the and -directions can be obtained as well
In the gas-kinetic scheme, the numerical fluxes at the Gauss quadrature point can be obtained by taking moments of the gas distribution function
where is provided by the integral solution of BGK equation Eq.(1) at the cell interface
and is the location of cell interface, is the particle velocity, is the trajectory of particles. For the second-order scheme, the second-order gas-kinetic solver for the three-dimensional flows GKS-Xu2 (); GKS-Xu1 () can be expressed as
where the coefficients in Eq.(2) can be determined by the spatial derivatives of macroscopic flow variables and the compatibility condition as follows
where the superscripts or subscripts of these coefficients are omitted for simplicity, more details about the determination of coefficient can be found in GKS-Xu2 (). With the second-order gas-kinetic solver Eq.(2), the second-order accuracy in time can be achieved by one step integration. In the one-stage gas evolution model, the third-order gas-kinetic solver has been developed as well GKS-high-1 (); GKS-high-2 (); GKS-high-3 (); GKS-high-4 ().
3 Fourth-order temporal discretization
In the further development of higher-order scheme, the formulation of one-stage gas-kinetic solver can become very complicated, such as the fourth-order method GKS-high-5 (), especially for multidimensional computations. Recently, a two-stage fourth-order time-accurate discretization was developed for Lax-Wendroff flow solvers, particularly applied for hyperbolic equations with the generalized Riemann problem (GRP) solver GRP-high () and gas-kinetic scheme GKS-high (). Such method provides a reliable framework to develop a three-dimensional fourth-order gas-kinetic scheme with a second-order flux function Eq.(2). Consider the following time-dependent equation
with the initial condition at , i.e.,
where is an operator for spatial derivative of flux. The time derivatives are obtained using the Cauchy-Kovalevskaya method,
Introducing an intermediate state at ,
the corresponding time derivatives are obtained as well for the intermediate stage state,
Then, the state w can be updated with the following formula,
It can be proved that for hyperbolic equations the above time stepping method Eq.(8) and Eq.(9) provides a fourth-order time accurate solution for at . The details of the analysis can be found in GRP-high ().
The semi-discretized form of finite volume scheme Eq.(2) can be applied for the above two-stage method . For the gas-kinetic scheme, the gas evolution is a relaxation process from kinetic to hydrodynamic scale through the exponential function, and the corresponding flux is a complicated function of time. In order to obtain the time derivatives of the flux function at and with the correct physics, the flux function should be approximated as a linear function of time within a time interval. According to the numerical fluxes at the Gauss quadrature points Eq.(5), the following notation is introduced
In the time interval , the flux is expanded as the following linear form
The coefficients and can be determined as follows,
By solving the linear system, we have
Similarly, for the intermediate state can be constructed. The corresponding fluxes in the - and -direction can be obtained as well. With these notations, the numerical scheme for the three-dimensional flows can be constructed, and the detail procedure of two-stage algorithm can be found in GKS-high ().
4 Spatial reconstruction
The above time evolution solution is based on the high-order initial reconstruction for macroscopic flow variables and WENO reconstruction WENO (); WENO-JS (); WENO-Z () is adopted for the spatial reconstruction. For the three dimensional computation, the reconstruction procedure for the cell interface is given as an example, and the stencil for reconstruction is given in Fig.1. The point value and and first-order derivatives at the Gauss quadrature points , need to be constructed. The detailed procedure is given as follows
According to one dimensional reconstruction, the cell averaged reconstructed value and cell averaged spatial derivatives can be constructed for the stencil shown in Fig.1.
With the one-dimensional WENO reconstruction in the horizontal direction, the averaged value and the averaged spatial derivatives and over the interval with can be given.
With one-dimensional WENO reconstruction in the vertical direction, the point value and spatial derivatives , and can be fully given at the Gaussian quadrature points .
In the computation, without special statement, the fifth-order WENO-JS reconstruction WENO-JS () is adopted for the flow with discontinuities and the linear scheme is used for the smooth flows to reduce the dissipation.
5 Numerical experiments
In this section, numerical tests for both inviscid and viscous flows will be presented to validate our numerical scheme. For the inviscid flow, the collision time takes
where and denote the pressures on the left and right sides of the cell interface. In the computation, and . For the viscous flow, we have
where is the viscous coefficient and is the pressure at the cell interface, and it will reduce to in the smooth flow regions. is the time step which is determined according to the CFL number, which takes in the computation. The reason for including artificial dissipation through the additional term in the particle collision time is to enlarge the kinetic scale physics in the discontinuous region for the construction of a numerical shock structure through the particle free transport and inadequate particle collision in order to keep the non-equilibrium property.
For the smooth flow, the WENO reconstruction can be used directly on the conservative flow variables. For the flow with strong discontinuity, the characteristic variables can be used in the reconstruction. Based on , where are the conservative variables, are the corresponding fluxes, and , the cell averaged and point conservative values can be projected into the characteristic field by , where is the matrix corresponding to right eigenvectors of . The reconstruction scheme is applied on the characteristic variables . With the reconstructed polynomials for characteristic variables, the conservative flow variables can be recovered by the inverse projection.
In the following, we present many numerical examples. Based on the same WENO-JS reconstruction, the results presented below have very high accuracy in comparison with other solutions from same reconstruction.
5.1 Three-dimensional Sod problem
This problem is a fully three-dimensional extension of the Sod problem. The computational domain is , and a sphere of radius separates two different constant states
The exact solution of this problem can be given by the one-dimensional system with geometric source terms
where the radial direction is denoted by , while is the radial velocity, is the number of space dimensions and . The two-stage fourth-order gas-kinetic scheme is used to solve this equation and the numerical results with cells are given as the reference solutions. In the computation, the uniform mesh with is used. The symmetric boundary condition is imposed on the plane with , and , and the non-reflection boundary condition is imposed on the plane with , and . The three-dimensional density and pressure distributions are given in Fig.3, and the density and pressure distribution along the line are given in Fig.3. The numerical results agree well with the reference solutions.
5.2 Flow impinging on sphere
In this case, the inviscid hypersonic flows impinging on a unit sphere are tested to validate robustness of the current scheme with different Mach numbers with . In the computation, a mesh shown in Fig.4 is used, which representing the domain in the spherical coordinate . The density distributions in the whole domain, Mach number and pressure distributions at the plane with , and pressure profiles along the surface of sphere with for Mach number and are shown in Fig.4, Fig.6 and Fig.6, where the shock is well captured by the current scheme and the carbuncle phenomenon does not appear Case-Pandolfi (). As comparison, the inviscid hypersonic flows impinging on a two-dimensional unit cylinder is also tested with mesh. The Mach number and pressure distributions for Mach number are shown in Fig.7.
5.3 Rayleigh-Taylor instability
The Rayleigh-Taylor instability occurs when an interface between two different fluids with different densities experiences a pressure gradient opposing the density gradient Case-rt (). This interface is unstable and any disturbance tends to grow, leading to the penetration of both fluids into each other. If the initial interface displacement is random, the Rayleigh-Taylor instability usually evolves into complicated turbulent mixing. In this case, the three-dimensional Rayleigh-Taylor instability is tested in a rectangular box with domain with a gravity pointed downward . Boundary conditions are applied at the four sides, while symmetric boundary conditions are applied at the top and bottom walls. The instability develops from the imposed single mode initial perturbation
where is the height of the interface and is the box width. The dimensionless parameters is the Atwood number
where are densities of heavy and light fluids, respectively. In the computation, the Atwood number is used. The the fluid interface from a single-mode perturbation at are given in Fig.8, and the density distribution at , and with are given in Fig.9. As expected, the heavy and light fluids penetrate into each other as time increases. The light fluid rises to form a bubble and the heavy fluid falls to generate a spike.
5.4 Lid-driven cavity flow
The lid-driven cavity problem is one of the most important benchmarks for numerical Navier-Stokes solvers. As shown in Fig.10, the fluid is bounded by a unit cubic and driven by a uniform translation of the top boundary. The monatomic gas with is used, such that there is no bulk viscosity involved. Early three-dimensional cavity-flow calculations were carried out by De Vahl Davis and Mallinson Case-Davis () and Goda Case-Goda (). In this case, the flow is simulated with Mach number and all the boundaries are isothermal and nonslip. Numerical simulations are conducted with three Reynolds numbers of and using meshes for the domain . The -velocity profiles along the vertical centerline line, -velocity profiles along the horizontal centerline in the plane with and the benchmark data Case-Albensoeder (); Case-Shu () are shown in Fig.11. The simulation results match well with the benchmark data.
5.5 Compressible homogeneous turbulence
The high-order gas-kinetic scheme is applied for the direct numerical simulations (DNS) of the compressible decaying homogeneous isotropic turbulence. The flow is computed within a square box defined as , and the periodic boundary conditions are used in all directions for all the flow variables DNS-1 (); DNS-2 (); DNS-3 (). In the computation, the domain is discretized with a uniform Cartesian mesh cells . A divergence-free random initial velocity field is generated for a given spectrum with a specified root mean square
where is a volume average over the whole computational domain. The specified spectrum for velocity is given by
where is the wave number, is the wave number at spectrum peaks, is a constant chosen to get a specified initial kinetic energy. The initial volume averaged turbulent kinetic energy and the initial large-eddy-turnover time is given by
The Taylor microscale Reynolds number and turbulence Mach number are given as
where is Taylor microscale
and can be determined from and with initialized and , and the dynamic viscosity is determined by
The time history of the kinetic energy, root-mean-square of density fluctuation and skewness factor for velocity slope are defined as
In the computation, , and , and the uniform meshes with , and cells are used. The iso-surfaces of criterion colored by velocity magnitude and the pressure distribution with at time are given in Fig.12. The time history of normalized kinetic energy , normalized root-mean-square of density fluctuation and skewness factor with respect to are given in Fig.13. The numerical results agree well with the reference data. With fixed initial , the cases with are tested, and the time histories of normalized kinetic energy are given in Fig.13 as well. With the increase of , the dynamic viscosity increases, and the kinetic energy gets dissipated more rapidly.
5.6 Taylor-Green Vortex
This problem is aimed at testing the performance of high-order methods on the direct numerical simulation of a three-dimensional periodic and transitional flow defined by a simple initial condition, i.e. the Taylor-Green vortex Case-Bull (); Case-Debonis (). With a uniform temperature field, the initial flow field is given by
The fluid is then a perfect gas with and the Prandtl number is . Numerical simulations are conducted with two Reynolds numbers and . The flow is computed within a periodic square box defined as . The characteristic convective time . In the computation, , and the Mach number takes , where is the sound speed.
The volume-averaged kinetic energy can be computed from the flow as it evolves in time, which is expressed as
where is the volume of the computational domain, and the dissipation rate of the kinetic energy is given by
The numerical results of the current scheme with mesh points for the normalized volume-averaged kinetic energy and dissipation rate with Reynolds numbers and are presented in Fig.15 and Fig.17, which agree well with the data in Case-Debonis (); Case-wang (). The iso-surfaces of criterions colored by velocity magnitude at and are shown in Fig.15 for and in Fig.17 for . The evolution of flow structure is evident, starting from large vortices and decaying into more complex structures. Different from many other higher-order methods, the current scheme has no internal degrees of freedom to be updated within each cell.
In this paper, based on the two-stage time stepping method, a fourth-order gas-kinetic scheme is proposed for the three-dimensional inviscid and viscous flow computations. With the three-dimensional WENO-JS reconstruction, a gas-kinetic scheme with higher-order spatial and temporal accuracy is developed. In comparison with the classical methods based on the first-order Riemann solver, for the same fourth-order accuracy in time the current scheme only uses two stages instead of four stages with the Runge-Kutta time-stepping technique. As a result, the two-stage GKS can be more efficient than the four-stage methods with the absence of two time consuming reconstructions. For the Navier-Stokes solutions, the current scheme doesn’t separate inviscid and viscous terms and they are treated uniformly from the same initial WENO-type reconstruction. The GKS can present very accurate viscous flow solutions due to its multidimensional flux function at a cell interface, where the gradients in both normal and tangential directions of flow variables participate in the gas evolution. The fourth-order GKS not only has the expected order of accuracy for the smooth flow, but also has favorable shock capturing property for the discontinuous solutions. Most importantly, the numerical tests clearly demonstrate that the current fourth-order scheme has the same robustness as the second-order one. The scheme has been tested from the smooth flows to the flows with discontinuities, and from the low speed to the hypersonic ones. For the three dimensional Navier-Stokes solutions, the current scheme is one of the state-of-art methods in the capturing of complicated flow structures.
The work of Pan is supported by the grants from NSFC (11701038) and China Postdoctoral Science Foundation (2016M600065), and the work of Xu is supported by Hong Kong research grant council (16206617, 16207715, 16211014) and the grants from NSFC (11772281, 91530319).
- (1) S. Albensoeder, H.C. Kuhlmann, Accurate three-dimensional lid-driven cavity flow, J. Comput. Phys. 206 (2005) 536-558.
- (2) M. Ben-Artzi, J. Falcovitz, A second-order Godunov-type scheme for compressible uid dynamics, J. Comput. Phys. 55 (1984), 1-32.
- (3) M. Ben-Artzi, J. Li, G. Warnecke, A direct Eulerian GRP scheme for compressible fluid flows, J. Comput. Phys. 218 (2006), 19-43.
- (4) M. Ben-Artzi, J. Li, Hyperbolic conservation laws: Riemann invariants and the generalized Riemann problem, Numerische Mathematik, 106 (2007), 369-425.
- (5) P.L. Bhatnagar, E.P. Gross, M. Krook, A Model for Collision Processes in Gases I: Small Amplitude Processes in Charged and Neutral One-Component Systems, Phys. Rev. 94 (1954) 511-525.
- (6) R. Borges, M. Carmona, B. Costa, W. S. Don, An improved weighted essentially non-oscillatory scheme for hyperbolic conservation laws, J. Comput. Phys. 227 (2008) 3191-3211.
- (7) J. R. Bull, A. Jameson, Simulation of the compressible Taylor-Green vortex using high-order flux reconstruction schemes, AIAA 2014-3210.
- (8) C. Cercignani, The Boltzmann Equation and its Applications, Springer-Verlag, (1988).
- (9) S. Chapman, T.G. Cowling, The Mathematical theory of Non-Uniform Gases, third edition, Cambridge University Press, (1990).
- (10) B. Cockburn, C. W. Shu, TVB Runge-Kutta local projection discontinuous Galerkin finite element method for conservation laws II: general framework, Mathematics of Computation, 52 (1989), 411-435.
- (11) B. Cockburn, C. W. Shu, The Runge-Kutta discontinuous Galerkin method for conservation laws V: multidimensional systems, J. Comput. Phys. 141 (1998), 199-224.
- (12) J.Debonis, Solutions of the Taylor-Green vortex problem using high-resolution explicit finite difference methods, AIAA Paper (2013) 2013-0382.
- (13) G. De Vahl Davis, G.D. Mallinson, An evaluation of upwind and central difference approximations by a study of recirculating flow, Comp. Fluids. 4 (1976) 29-43.
- (14) Z. Du, J. Li, A Hermite WENO reconstruction for fourth order temporal accurate schemes based on the GRP solver for hyperbolic conservation laws J. Comput. Phys. 355 (2018) 385-396
- (15) S. Gottlieb, C. W. Shu, Total variation diminishing runge-kutta schemes, Mathematics of computation, 67 (1998) 73-85.
- (16) K. Goda, A multistep technique with implicit difference schemes for calculating two- or three-dimensional cavity flows, J. Comput. Phys. 30 (1979) 76-95.
- (17) A. Harten, B. Engquist, S. Osher and S. R. Chakravarthy, Uniformly high order accurate essentially non-oscillatory schemes, III. J. Comput. Phys. 71 (1987) 231-303.
- (18) X.Y. He, R.Y Zhang, S.Y. Chen, G. D. Doolen, On the three-dimensional Rayleigh-Taylor instability, Phys Fluids 11 (1999) 1143.
- (19) J. Huang, K. Xu, P. Yu, A unified gas-kinetic scheme for continuum and rarefied flows II: multi-dimensional cases, Commun. Comput. Phys. 12 (2012) 662-690.
- (20) X. Ji, F.X. Zhao, W. Shyy, K. Xu, A Family of High-order Gas-kinetic Schemes and Its Comparison with Riemann Solver Based High-order Methods, J. Comput. Phys. 356 (2018), 150¨C173.
- (21) G.S. Jiang, C. W. Shu, Efficient implementation of Weighted ENO schemes, J. Comput. Phys. 126 (1996) 202-228.
- (22) J. Li, Z. Du, A two-stage fourth order time-accurate discretization for Lax-Wendroff type flow solvers I. hyperbolic conservation laws, SIAM J. Sci. Computing, 38 (2016) 3046-3069.
- (23) Q.B. Li, S. Fu, High-order accurate gas-kinetic scheme and turbulence simulation, SCIENTIA SINICA Physica, Mechanica Astronomica, 44 (2014) 278-284.
- (24) Q. Li, K. Xu, S. Fu, A high-order gas-kinetic Navier-Stokes flow solver, J. Comput. Phys. 229 (2010) 6715-6731.
- (25) W. Liao, Y. Peng, L.S. Luo, Gas-kinetic schemes for direct numerical simulations of compressible homogeneous turbulence. Phys Rev E, 80 (2009) 046702.
- (26) C. Liu, K. Xu, G.Z. Zhou, Cell Size Effect on Computational Fluid Dynamics: The Limitation Principle for Flow Simulation, arXiv:1706.08200v4 [physics.comp-ph].
- (27) J. Luo, K. Xu, A high-order multidimensional gas-kinetic scheme for hydrodynamic equations, SCIENCE CHINA Technological Sciences, 56 (2013) 2370-2384.
- (28) N. Liu, H.Z. Tang, A high-order accurate gas-kinetic scheme for one- and two-dimensional flow simulation, Commun. Comput. Phys. 15 (2014) 911-943.
- (29) X.D. Liu, S. Osher, T. Chan, Weighted essentially non-oscillatory schemes, J. Comput. Phys. 115 (1994) 200-212.
- (30) L. Pan, J. Li, K. Xu, A few benchmark test cases for higher-order Euler solvers, Numerical Mathematics: Theory, Methods and Applications, 10 (2017) 711-736.
- (31) L. Pan, K. Xu, A third-order compact gas-kinetic scheme on unstructured meshes for compressible Navier-Stokes solutions, Journal of Computational Physics, 318 (2016) 327-348.
- (32) L. Pan, K. Xu, Q. Li, J. Li, An efficient and accurate two-stage fourth-order gas-kinetic scheme for the Navier-Stokes equations, J. Comput. Phys. 326 (2016), 197-221.
- (33) L. Pan, K. Xu, A third-order gas-kinetic scheme for three-dimensional inviscid and viscous flow computations, Computers Fluids 119 (2015) 250-260.
- (34) M. Pandolfi, and D. D’Ambrosio, Numerical Instabilities in Upwind Methods: Analysis and Cures for the ”Carbuncle” Phenomenon, J. Comput. Phys. 166 (2001) 271-301.
- (35) W.H. Reed, T.R. Hill, Triangular mesh methods for the neutron transport equation, Technical Report LA-UR-73-479, 1973, Los Alamos Scientific Laboratory, Los Alamos.
- (36) R. Samtaney, D.I. Pullin, B. Kosovic, Direct numerical simulation of decaying compressible turbulence and shocklet statistics. Phys Fluids, 13 (2001) 1415-1430.
- (37) C. Shu, L. Wang, Y.T. Chew, Numerical computation of three-dimensional incompressible Navier.Stokes equations in primitive variable form by DQ method, Int. J. Numer. Meth. Fluids. 43 (2003) 345-368.
- (38) C. W. Shu, S. Osher, Efficient implementation of essentially non-oscillatory shock capturing schemes, J. Comput. Phys. 77 (1988) 439-471.
- (39) V.A. Titarev, E.F. Toro, ADER: arbitrary high order Godunov approach, J. Sci. Comput. 17 (2002) 609¨C618.
- (40) E. Toro, Riemann Solvers and Numerical Methods for Fluid Dynamics, Springer, (1997).
- (41) L. Wang, W. K. Anderson, Erwin, S. Kapadia, High-Order Discontinuous Galerkin Method for Computation of Turbulent Flows, AIAA Journal, 53 (2015) 1157-1171.
- (42) K. Xu, Direct modeling for computational fluid dynamics: construction and application of unfied gas kinetic schemes, World Scientific (2015).
- (43) K. Xu, A gas-kinetic BGK scheme for the Navier-Stokes equations and its connection with artificial dissipation and Godunov method, J. Comput. Phys. 171 (2001) 289-335.
- (44) K. Xu, Gas kinetic schemes for unsteady compressible flow simulations, Lecure Note Ser. 1998-03, Von Karman Institute for Fluid Dynamics Lecture (1998).
- (45) K. Xu, J. Huang, A unified gas-kinetic scheme for continuum and rarefied flows, J. Comput. Phys. 229 (2010) 7747-7764.
- (46) K. Xu, C. Liu, A paradigm for modeling and computation of gas dynamics, Phys. Fluids 29 026101 (2017).