Validating a Monotonically-Integrated Large Eddy Simulation Code for Subsonic Jet Acoustics Daniel Ingraham∗ and James Bridges† National Aeronautics and Space Administration, John H. Glenn Research Center, Cleveland, OH, 44135, USA The results of subsonic jet validation cases for the Naval Research Lab’s Jet Engine Noise REduction (JENRE) code are reported. Two set points from the Tanna matrix, set point 3 (M = 0.5, unheated) and set point 7 (M = 0.9, unheated) are attempted a a on three different meshes. After a brief discussion of the JENRE code and the meshes constructed for this work, the turbulent statistics for the axial velocity are presented and compared to experimental data, with favorable results. Preliminary simulations for set point 23 (M = 0.5, T /T = 1.764) on one of the meshes are also described. Finally, the a j ∞ proposed configuration for the farfield noise prediction with JENRE’s Ffowcs-Williams Hawking solver are detailed. Nomenclature x, y, z Cartesian coordinates, with origin at the nozzle exit and centerline, x indicating the axial coordinate, increasing downstream r Radial coordinate, r2 =y2+z2 ∆ Approximate mesh spacing L Length of the mesh region meant to resolve the turbulent nearfield D Exit diameter of the nozzle j U Fully-expanded jet exit velocity j U, U(cid:48) Mean and root-mean square of axial velocity, respectively T Fully-expanded jet exit temperature j T, T(cid:48) Mean and root-mean square of temperature, respectively T , T Startup and statistics-gathering simulation times, respectively 0 θ Polar angle, with θ =0 pointing upstream (i.e. the negative-x direction) y Estimate of viscous length scale in the nozzle boundary layer τ I. Introduction Considerable effort has been invested by the computational aeroacoustic community in the prediction of subsonic jet noise through large eddy simulation (LES)1–20 (see also the excellent review from Bodony21). Most of this work has been focused on using codes that employ high-order finite difference or finite volume schemes and structured grids in their discretization of the governing equations. A number examples of the application of unstructured LES codes to the problem of supersonic jet noise can also be found in the literature, e.g., the work done by researchers associated with Cascade Technologies22–26 and the Naval Research Laboratory (NRL).27–33 Comparably little work has been done to validate the effectiveness of unstructured LES codes at simulating the acoustics of subsonic jets. While not needed to capture the geometry of the simple nozzles considered in this paper, the flexibility provided by unstructured meshes allowsuserstousesmallcellsizesinlocalizedareasoftheflowwithoutaffectingthemeshdownstream,unlike structured meshes. Also, simulating the subsonic regime of jet flows approximating the takeoff condition is ofinteresttoNASA’sCommercialSupersonicTechnology(CST)project,withthegoalofprovidingacoustic predictionsthroughtraditionalFfowcs-WilliamsHawkingtechniques34andmoresophisticatedmethods(e.g., ∗AcousticsBranch,MS54-3,AIAAMember †AcousticsBranch,MS54-3,AIAAAssociateFellow 1of22 AmericanInstituteofAeronauticsandAstronautics the forthcoming Generalized Acoustic Analogy work of Leib et al.35). The combination of these needs has motivated addressing this gap in the literature. The Jet Engine Noise REduction (JENRE) application,36 developed by the NRL, is a computational aeroacoustics code designed to provide farfield acoustic predictions of jets using the LES approach to the simulation of turbulent flow. JENRE follows the Monotonically Integrated LES (MILES)37 philosophy of sub-grid scale modeling, relying on the inherent numerical dissipation of the code’s schemes to capture the effect of the unresolved turbulent scales will have on the resolved unsteady flow. Numerical dissipation is provided by blending low- and high-order finite element schemes through the Flux-Corrected Transport (FCT) approach. Proponents of this approach claim that the numerical dissipation effectively acts as an implicit model for the turbulent sub-grid scales, citing as evidence modified equation analysis38 and the results of applying such schemes to practical turbulent flows.39,40 Numerous finite element schemes are available in the version of the JENRE code used in this work – here, the second-order Taylor-Galerkin scheme41,42 has been chosen for spatial and temporal discretization, coupled with a lumped-mass Galerkin scheme through the (FCT) method.43 II. Test Case Parameters The flow parameters for the simulations described in this work are taken from the Tanna matrix,44 and areshowninTable1. ThenozzlegeometrychosenforthisworkwasSMC000,thebaselineroundnozzlefrom the Small Metal Chevron series developed at the NASA Glenn Research Center (GRC). Described in more detail by Bridges and Wernet,46 SMC000 is an axisymmetric nozzle with a nominal two-inch exit diameter that has been subjected to extensive experimental investigation at the GRC.45,46 Set Point 3 Set Point 7 Set Point 23 M 0.513 0.985 0.376 M =U /a 0.500 0.900 0.500 a j ∞ T /T 0.950 0.835 1.764 j ∞ p /p 1.197 1.860 1.102 0 j Re = ρjUjD 5.99·105 1.36·106 2.06·105 j µj Table 1. Set points used in this work. Two mesh generation programs were used to create the meshes used in this work: Pointwise47 and Gmsh.48 Thefirsttwomeshes,referredtohereasPW0andPW1,werecreatedwiththePointwiseprogram, while the third mesh, GM0, was created with Gmsh. All meshes used the same SMC000 nozzle geometry. The boundaries of the nearfield, sponge zone, and outer domain were also identical. Figure 1 shows a diagram of the relative locations of the nozzle geometry, nearfield mesh region, and sponge zone and outer mesh boundaries. Diagrams of PW0, PW1, and GM0’s mesh spacing distributions are shown in Figures 2, 3, and 4, respectively. Comparison of 2 and 3 shows that PW0 and PW1 differ in the element sizes in the region just downstream of the nozzle exit to about 12.5D , with PW0’s spacing in the former axial location j being twice as large as PW1’s, and increasing to the 0.04D value more quickly. Figure 5 shows a crinkle j slice of PW0’s z = 0 surface, where the surface of the slice has been constructed from the mesh element surfaces intersected by the z = 0 plane. Enlargements centered on the nearfield region for for this same crinkle slice are shown for the PW0 and PW1 meshes, respectively, in Figures 9 and 10. The element edges are colored in blue in these two figures to make the differences in mesh density more apparent. InspectionofFigure4showstheGM0meshusedapproximatelythesamemeshspacingsizesasthePW1 mesh; however,themeshspacinggrowsmoregraduallyfromtheminimum0.005D tothemaximum0.04D j j value. GM0 also includes a refined region around the nearfield with target mesh spacing of about 0.05D , j intended to propagate acoustic disturbances to a Ffowcs-Williams Hawking integration surface. The impact of this change can be seen by comparing Figure 6 to Figure 11. ThePW0andPW1meshescontainedapproximately15·106 and27·106 nodes,respectively,whileGM0, with its FW-H refinement region, has about 55·106. The nearfield region consists of a conical frustum of length L = 25D with a diameter of D at the nozzle exit and diverging at a rate of 4◦. The shape of the j j nearfieldwaschosentocapturethedevelopmentofthejet’sshearlayer,andwasinformedbytheshearlayer 2of22 AmericanInstituteofAeronauticsandAstronautics y 5D j sponge zone 5D D j j 40D j 4◦ x Dj 25D j nearfield region 25D j Figure 1. Diagram of the extent of the domain boundaries, including the size of the nearfield region and sponge zone. y ∆≈D j ∆(x=L)≈0.04D j ∆(x=L/16)≈0.01D j x ∆(x=L/4)≈0.04D j ∆(x=0)≈0.01Dj nearfield region length L=25D j Figure 2. Approximate mesh spacing distribution for the PW0 mesh. 3of22 AmericanInstituteofAeronauticsandAstronautics y ∆≈D j ∆(x=L/2)≈0.04D j ∆(x=L)≈0.04D j ∆(x=L/16)≈0.01D j x ∆(x=L/4)≈0.02D j ∆(x=0)≈0.005Dj nearfield region length L=25D j Figure 3. Approximate mesh spacing distribution for the PW1 mesh. y ∆≈D j ∆≈0.05D all along FW-H region j FW-H region ∆(x=L )≈0.04D nearfield j φ3Dj 11◦ x ∆(x=0)≈0.005D j 25D j Figure 4. Approximate mesh spacing distribution for the GM0 mesh. 4of22 AmericanInstituteofAeronauticsandAstronautics Figure 5. z=0 slice of the PW0 mesh. 5of22 AmericanInstituteofAeronauticsandAstronautics Figure 6. z=0 slice of the GM0 mesh, colored by instantaneous axial velocity. Figure 7. x=2Dj slice of the PW1 mesh, colored by instantaneous axial velocity. 6of22 AmericanInstituteofAeronauticsandAstronautics Figure 8. x=2Dj slice of the PW1 mesh, colored by instantaneous axial velocity, enlarged figure. Figure 9. Enlargement of the nearfield region of the z =0 mesh slice shown in Figure 2 for the PW0 mesh, with cell edges highlighted to make the mesh spacing distribution more clear. 7of22 AmericanInstituteofAeronauticsandAstronautics Figure 10. Enlargement of the nearfield region of the z=0 mesh slice shown in Figure 3 for the PW1 mesh, with cell edges highlighted to make the mesh spacing distribution more clear. thickness growth rates reported by Bridges & Wernet.46 Figures 6 and 11 seem to show that this choice was adequately large, as do Figures 7 and 8, x = 2D slices of the PW1 mesh. The nearfield length L is j comparableto15,16 orlongerthan4,5,7,17,19,20 thatusedbyotherresearchers, butthewidthof25D tan(4◦) j is considerably more narrow than similar simulations4,5,7,17,19,20 found in the literature (4D to 7.5D is j j typical). For all meshes, anisotropic tetrahedra were extruded off the interior nozzle walls to partially resolve the boundary layer, likely to be quasi-turbulent according to Bridges and Wernet.46 For PW0 and PW1, this was done using Pointwise’s T-Rex tool, while Gmsh’s transfinite interpolation was used for GM0. The spacing of these cells in the wall-normal direction is initially 5·10−4D , and grows at a rate of 1.2 until j the tetrahedra become isotropic. The length scale of the viscous sublayer for the flows studied here is approximately y =4.1·10−5·D (estimated from smooth flat plate boundary layer theory as described by τ j Schlichting49), and thus the simulations will not meet the “wall-resolved” requirement of ∆≈2y .50 τ The minimum mesh spacing in the nozzle nearfield region is 0.01D for the PW0 mesh, and 0.005D for j j the PW1 and GM0. Other subsonic jet LES examples in the literature use minimum mesh spacings ranging from about 0.002D 5,20 to 0.003D ,15,16 and past supersonic jet calculations done by NRL researchers used j j spacings as small as 0.0035D in shear layer.29,31 The mesh spacing increases gradually downstream to j 0.04D at the end of the nearfield region for both meshes. j Boundary conditions consisted of adiabatic no-slip for all nozzle walls, and the appropriate constant stagnation pressure and temperature for the set point under investigation for the nozzle inflow. JENRE’s characteristic-based farfield boundary condition was used for all other boundaries. As indicated in Figure 1, a sponge zone was used to damp any spurious reflections from these farfield. III. Results A. Set point 3 Forsetpoint3andthetwoPointwisemeshes,aconstanttimestepwaschosensuchthatthelocalCFLnumber would not exceed 0.65 anywhere during the calculation. This time step size was reevaluated periodically during the run (about every 100 steps). The GM0 runs used a constant timestep of 5.13 · 10−5D /U , j j correspondingtoaCFLofabout0.24. Eachsetpoint3runwasstartedfromauniformflowofzerovelocity and ambient pressure and temperature. Table 2 shows the amount of time the simulation was run before 8of22 AmericanInstituteofAeronauticsandAstronautics T ·U /D T ·U /D 0 j j j j PW0 215.0 300.0 PW1 747.0 547.0 GM0 202.0 145.0 Table 2. Set point 3 startup (T0) and statistics (T) times. gathering statistics (T ) and during the statistics-gathering process (T). The values for T ·U /D for the 0 j j PW0 and PW1 meshes comparable to the periods of 87.5D /U ,20 500D /U ,17 and 833D /U 19 found in j j j j j j the literature. GM0, the largest and thus most computationally expensive mesh, may benefit from more time. Figure 11. Instantaneous axial velocity non-dimensionalized by Uj for set point 3 run on the PW1 mesh. AfterperformingthesimulationswiththeJENREcodeasdescribedintheprevioussection,theresulting data was post-processed to obtain mean and RMS statistics of axial velocity, and then compared to Bridges andWernet’smeasurements.46 Figure11showstheinstantaneousaxialvelocityfromanarbitrarytimestep forsetpoint3onthePW1mesh,andseemstoindicatethatthewidthofthemesh’snearfieldregionislarge enough to capture the majority of the turbulent flow field. Figure12and13showthemeanandRMSofaxialvelocityalongthecenterlineofthejetasafunctionof the distance from the nozzle exit, and compare these values to the “consensus” dataset found in reference,46 which consists of a representation of an expected value for these quantities informed by experiments run by many researchers performed at multiple labs over a considerable length of time, including to a significant degree by Bridges and Wernet themselves. Figure 12 shows that little difference is observed between the results for the mean axial velocity U alongthecenterlineforthethreemeshes. Allsimulationspredictthelengthofthepotentialcore(theregion downstreamofthenozzleexitthatisunaffectedbytheshearlayer)accurately,andconsistentlyunderpredict U downstream slightly. More spread between the two simulations is seen in Figure 13’s root-mean square axial velocity U(cid:48) results. Both simulations match the experimental data reasonably well for x<10D , and j show significant dissipation for the region beyond about x>15D . Interestingly, the simulation performed j on the coarser PW0 mesh agrees better than the the two finer meshes. GM0’s RMS results show a bit more variation over x>10D , likely a result of the shorter statistics sampling time for this run. j 9of22 AmericanInstituteofAeronauticsandAstronautics 1.0 0.8 0.6 Uj / U 0.4 0.2 RCA/BridgesConsensus PW0,300TUj/Dj PW1,547TUj/Dj GM0,145TUj/Dj 0.0 0 5 10 15 20 25 x/D j Figure 12. Mean axial velocity U along centerline (r=0) for set point 3. RCA/BridgesConsensus PW0,300TUj/Dj 0.15 PW1,547TUj/Dj GM0,145TUj/Dj 0.10 Uj / (cid:48)u 0.05 0.00 0 5 10 15 20 25 x/D j Figure 13. Root-mean square of axial velocity U(cid:48) along centerline (r=0) for set point 3. 10of22 AmericanInstituteofAeronauticsandAstronautics