High-order accurate large-eddy simulations of compressible viscous flow in cylindrical coordinates

Highlights

• Development of a high-fidelity, parallel simulation tool for both DNS and LES.

• Accurate treatment of the centreline numerical singularity in compressible flows.

• Propose and validation of a new shock sensor based on a WENO smoothness indicator.

• Better performance in identification of the shock discontinuity.

• A 10% decrease of the computational cost using new proposed shock sensor.

Abstract

Dynamic interactions of shock waves and turbulent structures occur in a wide range of applications. The accurate simulation of turbulence in these flows requires a numerical scheme with minimal dissipation whereas the shock-capturing requires a dissipative scheme to stabilise the solution in the shock-wave discontinuous regions. These contradictory requirements make simulations of these flows extremely challenging. A compatible implementation of the solid boundary conditions at the point of incepting the acoustic waves and their internalisation into hydrodynamic instabilities is also challenging.

Here we present, a high-fidelity, massively parallel and accurate solver for both direct numerical simulations and large-eddy simulations of supersonic turbulent flows in cylindrical coordinates. A hybrid WENO/ high order central difference scheme is used for the spatial discretisation with a fourth-order five-stage 2N-storage Runge–Kutta for the time integration. The least square contraction of Lilly [1] is utilised for large-eddy subgrid-scale modelling. The solid surface boundary condition is implemented using the offset wall and ghost cell technique.

The numerical implementation is validated through a wide range of test cases from simple to more complex cases. The numerical results show good agreement with analytical solutions and available experimental results. In the large-eddy simulation of a supersonic under-expanded impinging jet, it is observed that the conventional Ducros sensor leads to a more dissipative hybrid scheme. In this paper, a new shock sensor based on a WENO smoothness indicator is proposed and it is demonstrated to perform better especially in the simulation of the supersonic under-expanded impinging jet by limiting the expensive WENO scheme to the discontinuous regions and reducing the computational cost by approximately 10%.

Introduction

In the last few decades, there have been significant developments in high-speed transport systems. In the operational condition of these systems, the Reynolds number is high (i.e. a wide range of turbulent scales) and Mach number is larger than unity (i.e. the existence of shocks is possible.), and to obtain an accurate and stable solution for the Navier–Stokes equations in such a complex flow is exceptionally challenging. To address this problem, different methodologies have been developed. One of the earliest approaches uses additional numerical dissipations to smooth the discontinuous regions [2], [3], [4], [5]. However, the numerical dissipation tends to filter out small-scale turbulence. This is undesirable for high-fidelity simulations where resolving all the dynamically significant turbulent structures is a necessary goal. Therefore, adding artificial numerical dissipations is considered to be undesirable in simulations of turbulent high-speed flows [6], [7].

There exist many numerical schemes with built-in numerical dissipation including but not limited to total variation diminishing (TVD) [8], flux-corrected transport (FCT) [9], essentially non-oscillatory (ENO) [10], weighted essentially non-oscillating (WENO) [11] and variations of these schemes [7], [12], [13]. The weighted essentially non-oscillating (WENO) scheme, developed originally by Jiang and Shu [11], has received considerable attention in the literature. The upwind discretisation in this method makes it robust and free of spurious oscillations in discontinuous regions¹; however, this method dissipates turbulent fluctuations. Different schemes based on the WENO have been developed over the last few years to address its shortcomings. Henrick et al. [12] observed that the fifth-order WENO method lost accuracy at the critical points (i.e. where the first derivative vanishes while the third derivative does not). They have proposed WENO5M where the first approximation of the stencil weights was calculated by using the weights in the original WENO scheme and then these weights were mapped to more accurate values, which leads to improved accuracy at a critical point. Borges et al. [13] also proposed the same framework where substantially larger weights were assigned to discontinuous stencils. Using this approach, a slightly better performance than WENO5M was achieved. More sophisticated developments of the WENO method have been proposed in recent years [7]. However, most of these methods are dissipative for turbulence simulations as they filter out small-scale turbulent structures while they capture the shock discontinuity accurately. Johnsen et al. [6] studied the capability of different approaches in the literature to model turbulent-shock phenomena and assessed them using a comprehensive range of test cases. It was shown that the WENO method with no additional modification is dissipative. However, hybrid WENO/central difference is not dissipative, captures the sharp shock profiles and resolves all length scales accurately. The only drawback of the hybrid scheme is the requirement for a precise definition of a shock sensor to activate either WENO or central difference schemes. Ducros et al. [16] proposed a shock sensor based on the fundamental differences of a shock and turbulenc. A shock discontinuity is associated with a substantial negative dilatation while turbulence is associated with strong vorticity. It should be noted that vorticity is generated downstream of the shock front as a result of a curved shock in the flow field, according to Crocco theorem [17], which also needs to be incorporated in this fundamental difference. Most of the previous studies only focused on the development of models for accurate capturing of shocks in laminar flows and in one and two-dimensional simple configurations [18]. The boundary conditions, especially the solid boundary conditions, are scarcely discussed in these studies. There are a few studies where these methods are integrated into a solver applicable to more practical configurations of shock/boundary-layer interaction [19], [20] and free/impinging jets [21], [22], [23], [24].

Therefore, one of the objectives of this study is addressing this issue with an efficient implementation of these methods into a high-accurate and massively parallel numerical solver. As discussed earlier, the hybrid approach relies on a sensor to determine in which regions the WENO scheme must be applied. There is still considerable controversy with respect to the definition of an accurate shock sensor. To minimise the computational cost of simulations, the expensive WENO method should be limited to the discontinuous regions while the stable solution is preserved. Therefore, another objective of this paper is to introduce a new shock sensor which is computationally cost optimised while resulting in a stable solution in the high resolution and parallel simulations of turbulent supersonic flows.

Another important aspect of high speed flows is the interactions of the flow with solid bodies. One of the consequences of this interaction is the receptivity process where acoustic disturbances are internalised into hydrodynamic instabilities. Receptivity is the key ingredient of the self-sustained oscillations observed in a plate with a cavity [25] and under-expanded supersonic jets [26]. This leads us to the third object of this study which is an accurate numerical implementation of the solid boundary conditions for the purpose of high-fidelity large-eddy simulation of supersonic impinging jets.

Having these goals in mind, this paper introduces an in-house parallel C++ code which is named ECNSS (Explicit Compressible Navier Stokes Solver). The code is developed for the solution of the compressible Navier–Stokes equations for turbulent flow where discontinuities exist due to shocks, interaction of the shock with turbulence is strong, and solid boundaries confine the flow.

The rest of the paper is organised as follows. The Navier–Stokes equations along with non-dimensionalised parameters are presented in Section 2. The large-eddy simulation subgrid model is discussed in Section 3. Sections 4 and 5 contain a brief overview of the Hybrid WENO/ Central difference and Runge–Kutta Schemes for the spatial and temporal discretisations, respectively. The inlet, outlet and adiabatic wall boundary conditions are explained in Section 6. Finally, the results for the test cases are presented in Section 7 which is followed by the conclusion in Section 8.