## Abstract

For jets, large eddy resolving simulations are compared for a range of numerical schemes with no subgrid scale (SGS) model and for a range of SGS models with the same scheme. There is little variation in results for the different SGS models, and it is shown that, for schemes which tend towards having dissipative elements, the SGS model can be abandoned, giving what can be termed numerical large eddy simulation (NLES). More complex geometries are investigated, including coaxial and chevron nozzle jets. A near-wall Reynolds-averaged Navier–Stokes (RANS) model is used to cover over streak-like structures that cannot be resolved. Compressor and turbine flows are also successfully computed using a similar NLES–RANS strategy. Upstream of the compressor leading edge, the RANS layer is helpful in preventing premature separation. Capturing the correct flow over the turbine is particularly challenging, but nonetheless the RANS layer is helpful. In relation to the SGS model, for the flows considered, evidence suggests issues such as inflow conditions, problem definition and transition are more influential.

## 1. Introduction

With recent dramatic increases in computing power, the use of large eddy simulation (LES) for practical engineering applications is becoming more feasible. As noted by DeBonis (2006), this has resulted in distinct areas of practical and rigorous LES. The former tends towards obtaining solutions on complex geometries at Reynolds numbers of practical interest. The approach tends to use robust, stable codes (often based on Reynolds-averaged Navier–Stokes (RANS) solvers) that use some form of upwinding (such as the Roe (1997) scheme) and hence naturally display dissipative qualities. In rigorous LES, more attention is paid to minimizing numerical errors and improved subgrid scale (SGS) modelling. Simulations are generally carried out at lower Reynolds numbers, giving a narrower range of turbulent scales and hence higher fidelity. Numerical methods for rigorous LES tend to be difficult to apply to real geometries. The drawback of a practical approach is that more empiricism is often needed to achieve a good solution, for example, in the level of numerical smoothing used.

For LES generally, as noted by Pope (2000), the numerical scheme can potentially have a strong influence on the results. As noted in the work of Ghosal (1996) and Kravchenko & Moin (1997), low-order numerical discretizations can be as influential as the SGS model. This means that it can be difficult to disentangle the effects of numerics and SGS modelling. To explore the effect of different numerics, one can omit the SGS model. Following Pope (2000), SGS model-free computations can be called numerical LESs (NLESs). If the code being used tends towards being dissipative, then omitting the SGS model can be helpful, otherwise excessive dissipation can occur. While the dissipative elements of the code are active and helpful in draining energy from the larger scales, it should not be supposed that they necessarily represent a particularly good SGS model. Also, for codes that do not exhibit numerical dissipation, SGS modelling can, if nothing else, be helpful for stability.

As implied by the work of Moulinec *et al*. (2005), the grid can have a strong solution influence. For example, meshes with lines connecting cell centres that are orthogonal to the cell faces can help reduce numerical diffusion and preserve kinetic energy (see Tucker *et al*. 2006).

The main pacing item for wall-bounded LES (Fureby 2004) is improved methods for modelling near-wall turbulence at high Reynolds numbers. Momentum transfer in the near-wall region is affected by streak-like structures. Pope estimates there to be 10^{8} streaks on a Boeing 777 wing at cruise. For high Reynolds numbers, it is too computationally demanding to resolve the streaks. Hence, it is necessary to model the effects of these near-wall structures and their transfer of momentum. A promising approach is the use of a RANS layer to cover over unresolved features and provide the mean flow characteristics of the wall boundary layer. Following Tucker & Davidson (2004), the RANS–LES interface can be specified at a location based on turbulence physics, *y*^{+}=60 is promising. This places the interface in the inner part of the logarithmic region, where the grid resolution demanding streaks and turbulence structures in the buffer layer are RANS modelled. Georgiadis *et al*. (2006) also successfully uses RANS–LES for a turbulent mixing layer, with the whole upstream boundary layer flow being RANS modelled.

To explore the issues discussed above, this work solves a jet flow for several numerical schemes using no SGS modelling. Also, several SGS models are investigated using the same scheme. Some wall-bounded flows are then investigated using a hybrid RANS–NLES approach.

## 2. Numerical methods

### (a) Numerical schemes

Table 1 shows six different numerical schemes considered. The top four codes, HYDRA (Lapworth 2004), FLUXp (Xia 2005), TBLOCK (Klostermeier 2008) and BOFFS (Jefferson-Loveday 2008) are similar to each other. They use second-order centred fluxes with an additional smoothing term containing a tunable parameter *ϵ*. The terms are shown in the right-hand column of table 1. HYDRA can be viewed as cell vertex based and unstructured. The flux at the control volume interface is based on the flux difference ideas of Roe (1997), where normally the inviscid flux at a control volume face, *Φ*_{f}, is expressed as(2.1)where *A*=∂*Φ*/∂*ϕ* (*ϕ* represents primitive variables) and the subscripts L and R represent variables interpolated to the control volume face based on information to the left- and right-hand side of the face. Following Morgan *et al*. (1991), *Φ*_{L} and *Φ*_{R} are simply taken as the adjacent nodal values, i.e. (*Φ*_{L}+*Φ*_{R})/2 represents a standard second-order central difference. The smoothing term is also approximated in a second-order fashion (shown in table 1), where and are undivided Laplacians evaluated at the node locations L and R. The term |*A*| involves differences (and summations) between the local convection velocity and speed of sound. As the Mach number (*Ma*) tends to unity, key terms in |*A*| become small. For the temporal discretization, an explicit five-stage Runge–Kutta scheme is used.

BOFFS is a multi-block, structured, curvilinear, second-order turbulence model test bed. The inviscid flux at the control volume face is expressed as in equation (2.1). The interpolation of the primitive variables (*ϕ*) is again shown in table 1. For the temporal discretization, a Crank–Nicolson scheme is used. FLUXp is a cell centred unstructured code using a modified Roe (1997) scheme. For the temporal discretization, the implicit second-order backward Euler scheme is used. TBLOCK is a structured cell vertex solver intended for turbomachinery use. The convective flux at the control volume face is the average of values stored at the four corner points of the surface. Hence, for an *i*=constant face,(2.2)The method is equivalent to the application of the two-dimensional trapezoidal rule, and it can be shown that the accuracy in space is second order. The temporal discretization is an explicit three-level scheme. The conserved variables are smoothed after each time step, where again *ϵ* is a user-defined value. The smoothing could be viewed with the filtering of Sagaut & Grohens (1999) in mind. CHYDRA (Ciardi *et al*. 2005) is a version of HYDRA that uses flow adaptive smoothing via additional terms to the right-hand side of (2.1). Where pointwise oscillations are detected, the smoothing is increased. VU40 (Tucker 2001) uses a one-legged Crank–Nicolson scheme with pure central differences and staggered grids. Unlike all the other codes noted above, VU40 is free of a second-order ‘smoother’ for shock capturing, again multiplied by tunable parameters. This term is likely to play an active role. VU40 is used to investigate the performance of different SGS models, which are described by Liu *et al*. (2008) and shown in table 2.

### (b) Ffowcs Williams–Hawkings surface

The Ffowcs Williams and Hawkings (FW–H) approach (Ffowcs Williams & Hawkings 1969) is used to calculate the farfield sound. As illustrated in figure 1, the control surface is placed around the nonlinear flow and noise sources. The pressure signal in the farfield can be found based on quantities on the control surface provided by the LES. Problems can occur due to the nonlinear flow crossing the downstream surface of the FW–H surface (called the closing disk problem). Figure 1 shows the parameter *L*, where 1*L* indicates that the entire closing disk is included and 0*L* indicates no closing disk.

### (c) Near-wall modelling

For the jet calculations using hybrid NLES–RANS, the modified HYDRA and FLUXp codes are used. HYDRA uses the *k*–*l* model of Wolfshtien (1969), whereas FLUXp uses the Spalart & Allmaras (1992) model. For both solvers, the RANS layer is used for *y*^{+}<60.0 (to save computational time, a developed RANS layer is generally frozen). For the RANS region, the exact wall distance, *d*, is required, whereas in the NLES region, the wall distance, *d*, is zero. Hence, a separate variable, , is defined, where in the RANS region, , and in the NLES region, . Following Tucker (2004), the near-wall RANS length scale involves Eikonal equation use (Sethian 1999). This is an exact differential equation for wall distance. Blending of the RANS and NLES regions is achieved by the addition of a Laplacian scaled by normal wall distance to the Eikonal equation. The Laplacian blended Eikonal can formally be classified as a Hamilton–Jacobi (HJ) equation. The capability to solve this equation has been added internally to both the HYDRA and FLUXp codes, but for the solutions presented here, these distributions are indirectly used. Figure 2 shows potential blending functions at the interface of the two regions.

For the compressor and turbine cases, nominally the RANS–NLES interface is set at *y*^{+}=approximately 200 at the leading edge of the blade. The RANS and NLES regions are again smoothly blended (see Klostermeier 2008). However, as noted by Tucker & Davidson (2004), the form of blending should not greatly influence solutions.

## 3. Sensitivity studies

### (a) Numerical scheme comparison and subgrid scale model comparison

As noted by Wagner *et al*. (2007), a recommendation of the use of any particular code for LES of acoustics cannot be made since a benchmark of codes has not yet been established. For SGS models, Spalart (2000) notes that the various models used by different researchers suggests a lack of maturity and/or indifferent progress in model development. Here, to consider the sensitivity of ‘LES’ to the numerical discretization, the six different structured and unstructured flow solvers in table 1 are compared. All the solvers use essentially the same 5×10^{6} cell mesh to simulate an isothermal, *Ma*=0.9 axisymmetric jet at a Reynolds number (*Re*) of 10 000. For the codes that display naturally dissipative qualities, no SGS model is used. Table 2 shows several LES SGS models used with the neutrally dissipative VU40 code. The derivative convention *u*_{i,j}=∂*u*_{i}/∂*x*_{j} is used. Linear and nonlinear models are compared using the same 1×10^{6} cell mesh and an isothermal plane jet case at *Re*=3000. Full details of the nonlinear models are given in Liu *et al*. (2008).

Figures 3 and 4 show results for the various numerical schemes and SGS models. Lines give predictions and the symbols give measurements. Figure 3 plots results for the different numerical schemes. Figure 3*a* shows the centreline velocity decay. Figure 3*b* shows the peak shear layer shear stress. The centreline velocity decay rates for the different numerical schemes are within the spread of the measurements. However, there is a difference in the potential core length of up to 20 per cent. For the peak shear stress, there is some variation in the axial location of the peak value. However by *x*/*D*>5.0, the decay rates are within 5 per cent for the different numerical schemes. The shear stress disparity in this area is also not just related to numerics. Inflow conditions and measurement uncertainties are factors. However, numerical diffusion will alter the effective Reynolds number, and this is most likely to have the strongest influence in the shear layers. Figure 4 plots results for the various SGS models. Again, figure 4*a* shows the centreline velocity decay and figure 4*b* shows the peak shear layer shear stress. For the various SGS models, it is difficult to distinguish any significant difference in the results for both sets of curves.

### (b) Acoustics

Figure 5*a* plots the overall sound pressure level (OASPL) against directivity at 100*D*_{j} from the nozzle. The HYDRA results are for the axisymmetric jet, the angle 0° is the downstream direction and 90° is the sideline direction. The symbols show the measurements of Tanna (1977). The lines show results for various extents of the FW–H closing disk, i.e. values of *L* between 0 and 1 (illustrated in figure 1). The variation of OASPL for different values of *L* is up to 8 dB. The trends for different *L* are similar, except for 0*L* (solid line), where there is a rise in OASPL for angles greater than 80°. An integration of the OASPL curves in figure 5*a* allows a single value to be plotted against *L*. This is shown in figure 5*b*, where the straight dashed line shows the integral of the measurements. The results show that the FW–H surface configuration can have a significant influence on the predicted sound level.

## 4. Coflowing geometry

Now a short cowl coaxial geometry more reminiscent of the nozzle found on a real jet engine is investigated using hybrid RANS–NLES with the modified HYDRA solver. The RANS layer is on the internal nozzle wall and cropped at the nozzle exit. The core and bypass streams are both isothermal, with jet exit velocities of 240 and 216 m s^{−1}, respectively. The geometry and flow are more complex than for the simple ‘jet’ just considered. There are two shear layers and contoured nozzle surfaces tapered to a single point. The mesh contains 12×10^{6} cells that are primarily concentrated in the near nozzle region. The experimental nozzle is small (*D*_{j}=0.018 m) and *Re* for the measurements and HYDRA computations are matched at 300 000 (based on the maximum jet exit velocity and bypass jet diameter). Figure 6 plots radial profiles of mean axial velocity (*a*), normal stress (*b*) and shear stress (*c*). Symbols show the measurements and the lines are the computations. The results are staggered to show the development of the jet at axial locations from *x*/*D*=1.0 to 10.0. Contour plots of the computations are shown in the background. The agreement between the measurements and the computations is encouraging, particularly for *x*/*D*>3.0. Nearer to the jet nozzle exit, there is some discrepancy in the Reynolds stress, particularly on the jet centreline, where the experimental normal stress is 60 per cent higher than the predictions. Figure 7 plots vorticity magnitude. Figure 7*a* shows qualitatively the outer shear layer growth, while figure 7*b* shows more of the downstream flow development to *x*/*D*=5.0.

## 5. Chevron geometry

Hybrid RANS–NLES is also used to investigate flow for the chevron nozzle shown in figure 8*a*. The RANS layer is on the internal nozzle wall and around the chevron features. The FLUXp solver is used with a 12×10^{6} cell mesh for an isothermal jet at *Ma*=0.9. Two chevron penetration angles are considered, 5.0° and 18.2°. The chevrons act to break up larger, more dominant structures that are responsible for downstream low-frequency noise. The break-up of these structures occurs through enhanced mixing, which is more violent for the high penetration angle than for the low angle. Figure 8*b*,*c* show instantaneous FLUXp ‘streamlines’ from the chevron tip, showing vortex roll-up a short distance from the tip. Figure 9*a* shows density gradient contours on the chevron tip and root planes. On the chevron root (9*a*(ii)), the contours spread out faster, whereas on the tip they are more constrained. The flow is, to some extent, forced out through the chevron root. Figure 9*b* plots the farfield sound at 40*D*_{j}. The lines show computations and the symbols show measurements. The solid line and filled symbols are for the 18.2° penetration angle. The dashed line and open symbols are for the 5.0° penetration. The increased chevron angle decreases the downstream noise but increases the sideline noise. The predictions show this in agreement with the measurements.

## 6. Compressor and turbine endwall flows

The endwall flows in both idealized compressor and turbine blade rows are investigated using hybrid NLES–RANS with the TBLOCK solver. The correct prediction of loss (see Denton 1993) is linked to capturing the development of the endwall boundary layer. Figure 10 shows coherent structures over the compressor blade and endwall, visualized by the *λ*_{2} criteria (see Jeong & Hussain 1995). The figure illustrates the general case set-up and interaction of the endwall flow with the blade suction surface. The RANS layer is included on the endwall surface only. Nominally, at the blade leading edge, *y*^{+}=200. Both the compressor and turbine exhibit secondary flow caused by the pressure gradient between the blades. Both cases have saddle points in front of the blade leading edges. This is the stagnation point on the endwall in front of the blade and results (in conjunction with the secondary flow) in vortex structures migrating across the blade passage. These vortices can ultimately lift off the endwall such that a thin new boundary layer develops behind the lift-off line. Also, in a compressor, the adverse pressure gradient tends to make the endwall boundary layer prone to separation. The complex flow physics requires a method that can handle large-scale vortical motions, boundary-layer separation, unsteady mixing processes and transitional events.

### (a) Compressor endwall flow

For the compressor flow, the *Re* (based on the chord and inflow velocity) is 230 000 with *Ma*=0.07. A 5×10^{6} cell mesh is used. Figure 11*a* plots the pitchwise mass-averaged spanwise variation of exit flow angle at 50 per cent chord downstream of the trailing edge. Here, the midspan is given by 0.5 and a span of 1 gives the endwall. The symbols show measurements, the full line the NLES and the dashed line the hybrid NLES–RANS results. Clearly, the RANS layer improves prediction of the exit flow angle. Figure 11*b* plots the spanwise variation of the total pressure loss coefficient at 50 per cent chord downstream of the trailing edge. The same line styles are used. The RANS–NLES gives uniformly better results than NLES. It captures the shape of the measurements with more accuracy and the absolute values more closely.

Figure 12 plots streamlines in the endwall boundary layer. Clearly, the NLES–RANS improves prediction of the saddle-point location. For NLES, the separation occurs too far upstream of the leading edge. The lack of turbulent fluctuations in the endwall boundary layer results in its inability to overcome the adverse pressure gradient and thus premature separation.

### (b) Turbine endwall flow

For the turbine endwall flow, *Re*=590 000 (based on chord and inflow velocity) and 5.4×10^{6} grid points are used. Figure 13*a* plots the total pressure loss coefficient against the axial chord. The symbols show the measurements, the solid line the NLES and the dashed line the hybrid RANS–NLES predictions. Again, the RANS layer is helpful. Figure 13*b* plots the spanwise variation of the total pressure loss coefficient at 9 per cent axial chord downstream of the trailing edge. The NLES results, shown by the full line, are in poor agreement with the measurements, represented by the symbols. The inclusion of the RANS layer helps, but the results are still fairly poor. This erroneous prediction hints at the inability of the hybrid calculation to capture certain passage flow phenomenon. This could be due to the interface between the RANS and NLES regions cutting through the developing vortex systems. In addition, large regions of the endwall can contain laminar flow. The change from turbulent inlet boundary layer to a laminar one cannot be accurately simulated by the RANS layer. While the RANS element is helpful in improving the results, care must be taken where there is complex flow physics that cannot necessarily be captured by hybrid methods.

Figure 14 plots streamlines in the endwall boundary layer. As with the compressor, unlike the hybrid simulation, the NLES does not capture the leading-edge saddle point and other features shown in the frame (*a*) experimental observations. It should be noted that, for these blade flows, simple RANS modelling can give close agreement with measurements.

## 7. Conclusion

Increasing computer power is enabling LES to be applied to industrial flows beyond the simpler geometry canonical flows considered in the past. The techniques used when solving the former can be quite different to the latter, drawing a distinction between ‘practical’ and ‘rigorous’ LES. Here, a practical LES-type approach is used, which tends to focus on more complex geometries at higher Reynolds numbers using more robust codes. Several such solvers were used with no SGS model to investigate jet flows. Following this, several SGS models were tested using one numerical scheme, also for a jet flow. The results showed that the choice of SGS model was relatively unimportant and that, for schemes which display dissipative qualities, SGS model omission can yield useful results.

Following this, solutions with no SGS model but a RANS layer were made to investigate complex nozzle geometries and compressor and turbine endwall flows. For the cases considered it was, in a sense, demonstrated that, for practical flows, inlet conditions and problem definition are more important than the SGS model. The use of a near-wall RANS layer was shown to help give better agreement with measurements. However, for complex flows, such as the turbine endwall, where for example relaminarization can occur, the precise effects of using a near-wall RANS layer are unclear. Hence, while the use of the RANS layer is clearly helpful, some complex flow physics cannot be captured with a hybrid method and user judgement is needed.

## Acknowledgements

The authors would like to thank John Denton, Richard Jefferson-Loveday and Rolls-Royce Plc for the use of their codes (TBLOCK, BOFFS and HYDRA, respectively). P. G. T. would like to gratefully acknowledge funding under the Royal Society Industry Fellowship scheme.

## Footnotes

One contribution of 16 to a Discussion Meeting Issue ‘Applied large eddy simulation’.

- © 2009 The Royal Society