A Parametric Study on the LES Numerical Setup to Investigate Fan/OGV Broadband Noise

: In the present paper, large eddy simulations are performed to study two different mechanisms of Fan/OGV broadband noise: airfoil self-noise and turbulence interaction noise. Firstly, the current study focuses on the prediction of airfoil self-noise from a thin plate with a sharp trailing edge and a chord-based Reynolds number of the order of 10 6 . The boundary layer is tripped to trigger transition to turbulence, and a parameter study is performed to study the inﬂuence of the near-wall modeling, grid topology and reﬁnement in the near-wall and wake regions, the spanwise domain extent, and the tripping method. Empirical and analytical models, as well as available DNS data are used for validation purposes. Secondly, the interaction noise from a thin plate impinged by an incoming synthetic turbulent ﬂow is studied. For both cases, far-ﬁeld acoustic spectra are compared to Amiet’s models for leading and trailing edge noise showing a good agreement.


Introduction
The increase in airports capacities and stringent noise regulations have turned noise pollution in the vicinity of airports into a significant challenge for aircraft manufacturers. In the last decades, jet noise has been significantly reduced by increasing the bypass ratio of turbofan engines. Additionally, the development of passive noise control technologies, such as acoustic liners and the use of convenient blade/vane counts have led to significant tonal noise reduction. Recently, modern ultra-high bypass ratio (UHBR) engines have allowed for further noise reduction while maintaining the fuel consumption reduction trend. Fan broadband noise becomes a major contributor to the total radiated noise. Thus, the understanding of the broadband noise has become a relevant research topic. Fan broadband from an aero-engine is generated by a number of mechanisms, such as turbulence ingestion, tip leakage vortex, rotor/stator self-noise, and turbulent fan wake impingement onto the outlet guide vanes (OGV). The last two mechanisms, also known as trailing edge noise (TEN) and turbulence interaction noise (TIN, or leading edge noise), respectively, illustrated in Figure 1, are the object of the present numerical study.
Rotor broadband self-noise is generated by the scattering of turbulent kinetic energy from the boundary layer into acoustic waves as the flow encounters the singularity at the airfoil trailing edge. The estimation of airfoil self-noise requires an accurate prediction of the turbulent boundary layer development close to the trailing edge. For some airfoil configurations, at low to moderate Reynolds numbers, boundary layers may still be transitional as they reach the trailing edge. This is usually the result of large laminar flow regions around the airfoil, leading to self-sustained oscillations and Tollmien-Schlichting instability waves developing on the suction or pressure side of the airfoil. These waves may then interact with the trailing edge and generate additional noise. In order to study trailing edge noise at higher Reynolds numbers, such instabilities should be avoided. A possible solution to force the transition to turbulence is to trip the boundary layer [1,2]. The TIN results from a direct impingement of the rotor turbulent wakes on the stator leading edge. Turbulence coming from the rotor induces load fluctuations on the vanes, which produces broadband noise.
Different numerical approaches can be used to investigate airfoil broadband noise [3,4]. Direct numerical simulation (DNS) is convenient, as no turbulence modeling is required. However, it is still limited to simple cases due to its prohibitive CPU cost. The reduction in the computational cost requires the use of some level of turbulence modeling. In this work, large eddy simulation (LES) is used, in which only the smallest and poorly energetic turbulent eddies need to be modeled. The acoustic propagation can be directly resolved by LES. Alternatively, turbulent flow parameters can be extracted from LES to feed analytical broadband noise models.
The main objective of this paper is to investigate the numerical requirements for a correct description of the TEN and TIN mechanisms with LES. Comparisons to exact analytic models are preferred to avoid any experimental or numerical uncertainty in the reference data. Consequently, each noise mechanism is investigated separately on a dedicated flat-plate configuration. Flow conditions are chosen to be representative of the fan stage aero-engine at approach conditions with a significant chord-based Reynolds number (Re = 10 6 ) and a freestream Mach number of about M = 0.3. Beyond the investigation of the numerical parameters, this study is also the opportunity to clarify the dependency of broadband noise on the main physical parameters.

Analytical Models
The analytical models used in this study are based on Amiet's work for single airfoils. The model assimilates thin airfoils, such as turbofan rotor blades or stator vanes, to flat plates with no thickness or camber and zero angle of attack.
Amiet's TEN model [5] describes the scattering of pressure fluctuations from an incident boundary layer at the sharp trailing edge [6]. Assuming a large span to chord aspect ratio, the power spectral density (PSD) of the acoustic far-field, for an angular frequency ω and an observer at (x 1 ,x 2 ,x 3 ), is given by where k is the acoustic wave-number, k = kc/2, c is the chord length, L is the span length, k 2 is the spanwise wave-number, S 0 represents a convection-corrected distance, and u c is the convective velocity. Π 0 is a statistical function related to the wall-pressure spectrum φ pp slightly upstream of the trailing-edge, assuming homogeneous boundary layer turbulence. l z is the associated spanwise correlation length where γ 2 is the spanwise coherence between points separated by a distance η 2 . In Equation (1), I is the aeroacoustic transfer function that can be considered as the sum of the contribution of the main scattering from the trailing edge I 1 and the leading edge scattering correction I 2 . Detailed expressions can be found in [6]. The TIN model considered here was first proposed by [7]. The incident wake impinging on the leading edge is assumed to be frozen turbulence and is described as a sum of spatial Fourier modes in the streamwise and spanwise directions. Assuming a large span to chord aspect-ratio, the resulting far-field sound at an observer position in the midspan plane (x 1 ,0,x 3 ) is given by where k 1 is the streamwise wave-number, Φ ww is the PSD of the upwash velocity fluctuations, and L(x, k 1 , k 2 ) is the total aeroacoustic transfer function. For a frozen turbulent velocity field k 1 = ω/u 0 .

Trailing Edge Noise (TEN)
The simulations address the flow around a flat plate with a sharp trailing edge at a free-stream velocity u 0 = 100 m/s, in atmospheric conditions (ρ ∞ = 1.177 kg/m 3 , T ∞ = 300 K, µ ∞ = 1.81 × 10 −5 Pa.s). The chord is c = 0.1 m and the thickness is equal to c/100. The Reynolds number based on the chord is Re c = 10 6 . The computational domain size is 20c in the x-direction and 25c in the y-direction. The spanwise extension is 10% of the chord length, which represents about five times the boundary layer's thickness at the trailing edge. An additional simulation was performed with a span length equal to the chord length in order to assess the effect of confinement on the results. Unless mentioned otherwise, periodic boundary conditions are used in the spanwise direction. Non-reflective Navier-Stokes characteristic boundary conditions (NSCBC) [8] are used at the inlet and outlet planes of the domain, as well as at the upper and lower boundaries. The AVBP explicit unstructured compressible LES solver developed at CERFACS [9] is used to solve the governing equations, with the SIGMA sub-grid scale model [10], and the two steps Taylor Galerkin TTGC scheme [11], which is a third-order convective numerical scheme in space and time.
A parameter study is performed in order to study the influence of the near-wall modeling, the grid topology and refinement in the near-wall and wake regions, the spanwise domain extent and boundary conditions, and the tripping methodology. The nomenclature for the different cases is as follows: [ For all the cases, the mesh is mostly composed of unstructured tetrahedral elements and is carefully clustered near the airfoil surface to achieve an adequate spatial resolution. Near the wall, two types of elements are compared : tetrahedral and prismatic. For the latter, 8 layers of prisms are superimposed. The smallest wall normal spacing corresponds to the first layer near the wall and increases progressively with a geometric ratio of 1:1. In all the cases, a square trip is used and the mesh around the geometric trip is tetrahedral in order to avoid very small cells, which would lead to inconveniently small timesteps.

Influence of the Mesh Type and Refinement at the Wall
The influence of the mesh type near the wall is studied by comparing tetrahedral and prismatic meshes. The trip height corresponding to H4 is used for this comparison since it gives the best transition in WM cases (as it will be shown in the next section). Figure 2 (left) shows u + , the mean streamwise velocity u normalized by the shear velocity u τ , as a function of the distance to the wall y + , in wall units, at x/c = 70%. LES results are compared to the analytical solution of VanDriest [12]. The WR case (GeomX10H4-P01W30) shows the best agreement as expected, with a slightly larger value of u + in the buffer region when compared to the WM cases. In the WM cases, for the same value of y + = 25 at the wall and the same CPU cost, a better agreement is found for the prismatic mesh when compared to the tetrahedral one. In order to reach a similar agreement, the tetrahedral mesh needs to be 2.5 times finer at the wall, which leads to y + = 10. The influence of the mesh near the wall on the wall pressure spectrum (WPS) in the vicinity of the trailing edge is also analyzed. Figure 3 shows the WPS obtained from the LES data, for prismatic and tetrahedral cases, compared to the empirical solution proposed by Goody [13]. The numerical WPS are normalized by u 0 /τ 2 w δ and plotted against the reduced frequency ωδ/u 0 . It should be noticed that the WPS is quite sensitive to the mesh type and refinement close to the wall. A very good agreement with the empirical model is obtained in the WR case (GeomX10H4-P01W30), while some discrepancies can be found for the WM cases. For the tetrahedral case (GeomX10H4-T25W30), more discrepancies are visible, particularly at high frequencies.

Influence of the Tripping Methodology
Two tripping methods are compared in this study. The first one corresponds to a geometric trip on the plate surface, with a rectangular step shape, which requires a dedicated mesh. The second one involves a source term introduced in the streamwise momentum and energy equations, which has been designed to mimic the presence of a roughness element in the flow [14]. The root-mean-square (RMS) velocity profiles u rms for the two tripping methods are compared in Figure 2 (right) for WR simulations at a streamwise position corresponding to Re θ = 1968 (Reynolds number based on the momentum thickness of the boundary layer θ). They are also compared to the DNS data from [15]. Although the large geometric trip ("H4") overestimates the turbulent intensity in the boundary layer, a very good agreement with the DNS data in terms of shape and amplitude is obtained using the source term. Thus, the source term allows for a correct transition without over-estimating the turbulence intensity.

Influence of the Trip Height
The influence of the trip height is studied in WM and WR cases. Three heights are considered for each case H trip δ * = 1, 2, and 4. The instantaneous flow topology is shown in Figure 4 for the trip heights H1 (a) and H4 (b) using iso-surfaces of Q-criterion from the LES computations. A transitional and turbulent boundary layer is found on both sides of the plate. The transition is triggered by the geometric trip with small structures appearing close to the reattachment point. Although similar turbulent structures are observed for both trip heights in the WR cases (not shown here), the WM simulation with the larger trip (H4) shows finer structures than that with the H1 trip. Note that for the WM cases, the small trip height (H1) is of the order of the first layer of prisms. The development of the turbulent boundary layer is then analyzed in Figure 5 by comparing the normalized velocity profiles and the WPS near the trailing edge for the different trip heights. Identical velocity profiles are obtained for the WR cases with the three trip heights, in very good agreement with the analytical solution. However, as it can be seen in Figure 5 (left) for the WM cases, the simulation corresponding to the largest trip "H4" shows a better agreement with the analytical solution. Thus, for WM simulations, the trip height must be sufficiently larger than the boundary layer displacement thickness in order to ensure an adequate boundary layer transition. This effect on the boundary layer transition can also be found by comparing the WPS profiles in Figure 5 (right). Only slight differences are obtained between the WR cases. However, for the WM cases, the highest trip (H4) shows a reduction in the WPS amplitude by about 10 dB/Hz compared to H1 and H2 trips at low to moderate frequencies, but an increased intensity at high frequencies, and presents the closest spectrum to the ones obtained in WR cases.

Influence of the Spanwise Extent and Boundary Condition
The influence of the spanwise extent and boundary condition on the spanwise correlation length l z is examined in Figure 6 (right) through a set of WM simulations corresponding to GeomX10H4-P25W30. Periodic ("PERIO") and symmetry ("SYM") boundary conditions are considered and two span extents (10% and 100% of the chord length) are investigated. In addition, l z obtained from a WR case is given in Figure 6 (left) as a reference. The LES results are compared to Corcos [16], Guedel [17], Efimtsov [18] and Salze [19] empirical models. Except for the models of Efimtsov [18] and Salze [19] that give slightly close results for l z , large differences can be found with the other models. For the WR case, a good agreement with the model of Salze [19] is observed for high frequencies and a slight over-prediction at low to mid frequencies is found due to the domain size. The model of Salze [19] can, thus, be used as a reference for the WM cases. Figure 6 (right) shows that periodic boundary conditions lead to a relatively better agreement with the model of Salze [19] at high frequencies ( f > 4 kHz) compared to symmetric boundary conditions, while a slightly higher over-prediction at low to mid frequencies can be noticed. When periodic boundary conditions are used, a perfect correlation is imposed between the boundaries, but the normal velocity fluctuations of turbulence are allowed. This explains why l z is overestimated at low frequencies, which correspond to length scales of the order of the spanwise extent of the domain, and are correctly estimated at intermediate and high frequencies, which are mainly controlled by turbulence. The comparison between the two span lengths is also showed in Figure 6 (right). Even if the spanwise extent of the computational domain (L LES = 10 mm) appears sufficient when compared to the maximum correlation length (l z,max = 3.5 mm), a larger extent improves the correlation length at lower frequencies. It can be observed that the spanwise correlation length l z is significantly reduced, particularly at low and mid frequencies, when the large span is considered. The range where l z is over-predicted is restricted to lower frequencies.

Influence of the Mesh in the Wake Region
For a fan/OGV stage, the wake formed behind the fan reaches the leading edge of the stator and generates rotor-stator interaction noise, which is the dominant source in such configurations. Thus, the numerical setup should allow for correct modeling and propagation of the fan wake in the inter-stage. After some distance downstream of the trailing edge of the fan blade, the wake reaches an asymptotic state. Once this length scale is known, the wake can be decomposed into three different regions: "near wake", "intermediate wake", and "far wake". The asymptotic far wake stage is reached when all turbulent flow statistics follow universal distributions that are independent of the trailing edge conditions. Here, numerical results obtained from the LES cases are compared to empirical models or universal distributions in each wake region.
Three mesh refinements are considered in the wake region, denoted "W100", "W30", and "W15". The WM case corresponding to the mesh type and refinement at the wall "P25" is used for all three cases. Figure 7 (left) shows the streamwise development of the center-line velocity (U c ) for the three mesh refinements in the near wake region. The shear velocity u τ0 is computed at the trailing edge of the plate and the kinematic viscosity ν is equal to its averaged value at each location. The LES results are compared with the near-wake analytic model proposed by [20]: where y * = u τ0 y ν , ψ = y/g(x), E 1 (ψ) = ∞ ψ e −t t dt, the function g is defined as g(α)[ln(g(α)) − 1] = κ 2 α, and γ = 0.5772157 is the Euler-Mascheroni constant. κ and B are model constants equal to 0.418 and 5.5, respectively. It can be noted that beyond a certain streamwise position (x * = 2 × 10 3 ), a good agreement is found between the LES results and the analytical solution, although U c is under-predicted for the coarsest mesh "W100".
The intermediate wake characteristics are analyzed for the "W30" case by comparing the Reynolds stress profiles for different streamwise positions in Figure 7 (right).
dy is the momentum thickness in the far-wake, with U the longitudinal mean velocity, and u 0 the constant freestream velocity outside the wake. In this region, the flow is expected to follow the asymptotic behavior given by [20] where η = y/b and U d is the velocity deficit. The profiles are plotted against the normal coordinate y normalized by the local wake half-width b. For η < 1, the normalized profiles are quite similar and fit well with the analytical solution. For larger values of η (i.e., y > b), the analytical solution predicts a turbulent intensity that tends to zero, while the Reynolds stress from the numerical results approaches its value outside of the wake region. The far-wake properties are investigated in Figure 8. The asymptotic behavior is reached approximately beyond x/θ = 350. The results in the far-wake can be compared to the half-power laws for the half-width b and maximum velocity deficit U d , which are given by In the case of the refined meshes "W30" and "W15", the power laws for b and U d are well predicted and the product U d b/u 0 θ approaches the predicted constant value 0.9394 [20]. For these cases, ν t /u 0 θ approximately tends to 0.032 and G max is slightly less than 0.0487 (see [20]). The mesh resolution "W100" seems to be too coarse in the wake region. Thus, it can be concluded that the mesh resolution "W30" is sufficient to properly predict the flow in the wake region beyond a certain streamwise position, with an acceptable computational cost.

Acoustic Far-Field Predictions
The acoustic far-field predicted using the LES results from the GeomX10H4-P25W30 and GeomX10H4-P01W30 simulations (both providing 13 points per wavelength in the far-field and a total number of points of 3.1 and 13 millions, respectively) is compared to Amiet's TEN model (see Equation (1)). For the analytical predictions, the convection speed u c is computed as the slope of the phase of the wall-pressure cross-spectrum between two streamwise locations in the vicinity of the trailing edge. The calculated values are between u c /u 0 = 0.63 and u c /u 0 = 0.7. The wall-pressure PSD in the vicinity of the trailing edge φ pp is extracted from the LES computation and is fitted to a Goody model, which can be used as an input in the TEN model. The spanwise integral lengthscale, l z , is modeled using the analytical models shown in Figure 6. For the direct numerical prediction, both the turbulent flow and the acoustic wave propagation are modeled by the LES solver. A comparison between the LES predictions and the analytical solutions of the PSD of the pressure fluctuations in the far-field, at 90°and a distance of 8c from the trailing edge of the plate, is shown in Figure 9. The mesh cut-off frequency can be observed at around 15 kHz, the precise value depends on the mesh resolution. The analytical prediction is very sensitive to the model used for l z . Both the WR and WM results show a good agreement with the model of Salze [19] at mid-to-high frequencies, but seem to overestimate the low frequency levels. This can be explained by the domain confinement in the spanwise direction, as shown for l z in Figure 6.

Turbulence Interaction Noise (TIN)
Similar computational domain and numerical setup to the TEN simulations are adopted here. The only difference is related to the inlet boundary condition, where synthetic turbulence is injected, based on a summation of Fourier modes [21]. The mean axial flow speed is u 0 = 100 m/s. At a distance of 0.05c upstream the leading edge, the turbulent spectrum can be fitted to a Von-Karman spectrum with a turbulence intensity of 4% and an integral length scale of 5 mm (5% of the chord length). The cut-off wavelength is chosen as twice the largest cell size according to Shannon's principle. The number of Fourier modes that compose the random velocity perturbation is 50. Turbulent kinetic energy follows the Passot Pouquet spectrum [22].
A total of four different grids are considered with 8 layers of prisms imposed on the plate surface. The mesh refinement at the wall is the same for all cases, and provide a consistent description of the RMS pressure along the plate (not shown). The only difference consists in the far-field cell size normalized by the chord length, equal to 1.7%, 2.1%, 3.6%, and 5.5% for MESH4, MESH3, MESH2, and MESH1, respectively. Figure 10 shows the numerical power level spectra (PWL) for the different meshes at a distance of 8c (left), and the PSD of pressure fluctuations for MESH4 (right) for a probe at 90°above the plate leading edge. For comparison with Amiet's TIN model, the PSD Φ ww and the spanwise correlation length l z (ω) of the normal velocity fluctuations are extracted from the LES simulation slightly upstream of the leading edge and used as inputs for the model. Apart from slight differences at low frequencies, a very good agreement can be observed with the theory for all cases up to 4 kHz. All these meshes have the same structure and the same minimum cell size near the wall where the noise is generated. The far-field grid refinement influences the propagation and consequently the cut-off frequency of the spectrum, such that the most refined mesh can propagate a larger range of frequencies. One can infer from Figure 10, for the specific numerical setup, that 13 points per wavelength are required to accurately capture acoustic waves up to a frequency of 20 kHz, which represents the audible limit. This valuable information will be used in the future to set up aeroacoustic-oriented LES simulations. Finally, it can be observed in Figure 10 (right) that the numerical pressure PSD nicely reproduces the lobes of the analytical solution, which are due to the non-compact acoustic sources on the plate surface.

Conclusions
A parametric study on the effects of several numerical parameters on the flow development around a flat plate and the associated noise generation was conducted using LES simulations. The main noise mechanisms were investigated separately: trailing edge noise and turbulence interaction noise. The upstream flow properties and the geometric parameters of the plate were chosen to approach as close as possible the assumptions made in Amiet's acoustic models. Comparison with analytical models and DNS data show very good results for the boundary layer characteristics and the far-field noise prediction. These results allow the definition of an adequate numerical set-up in order to prepare LES simulations of the broadband noise generated by an UHBR turbofan fan stage.

Conflicts of Interest:
The authors declare no conflicts of interest.

Abbreviations
The following abbreviations are used in this manuscript: