Spectral/hp Methodology Study for iLES-SVV on an Ahmed Body

This work focuses on the correlation study between a computational and physical model of an Ahmed Body with slant angle of 25 ◦ , which generates a complex ﬂow behaviour over the slant and back, with two vortices being generated from the side combined with separation on the slant. Physical results are from a wind tunnel test, performed by Strachan et al. [12] considering moving ground and Reynolds number of 1.7M, based on the length of the body. CFD simulations were performed using the code Nektar++, which is an open source, spectral/hp element high-order solver, which methodology combine both mesh reﬁnement (h), with higher polynomial order (p) for higher ﬁdelity modelling. It employs an implicit type turbulence model using a Spectral Vanish Viscosity (iLES-SVV) model, which works as a ﬁlter for high frequencies. Same physical test conditions and tunnel test section were also considered, over a total time of 4 convective lengths, with same Reynolds number of 1.7 Million from reference experiments. Considering the was a maximum not consider the upper Comparing the Spectral/hp element LES-SVV case from literature, the agreement with the experimental drag coefﬁcient has been improved, reducing the gap from 45 to 16%. For the lift coefﬁcient the maximum difference between the simulation results compared to experimental data is only 3%. There is also a good agreement between the LDA measurements on the end of the body with the results from the simulation. It is possible to observe a more intense vortex core on the simulation results, as compared to experimental data, which might well be explained by the upper support used to ﬁx the Ahmed body in experimental test, which weakens the vortices. The methodology shows promising results against the open literature once an appropriate validation study has been undertaken. Despite the relatively coarse resolution adopted the results are encouraging. Having identiﬁed an appropriate resolution, we will next consider other slant angles, to see how well these correlate with the experimental studies.


Introduction
Among all automotive bluff bodies in literature, the most studied one is the Ahmed Body. It was first proposed by Ahmed et al. [2], based on previous work from Morel [7], which was the first to study the behavior of slanted bluff bodies. The Ahmed body was designed to have shape similar to road vehicles and generate their main flow features, such as stagnation and separation points. The main dimensions of the Ahmed body are highlighted on Fig. 1.
Based on the results found by Ahmed et al. [2] on the variation of the slant inclination angle, Huminic and Huminic [4] states that three different flow configurations are found: from 0 to 12.5 • , the airflow over the angled surface remains fully attached before separating from the model when it reaches the vertical surface of the back end. The flow from the angled section and the side walls produces a pair of counter rotating vortices, which continue downstream; from 12.5 to 30 • , the flow over the angled section becomes highly complex. Two increased counterrotating lateral vortices are shed from the sides of the angled section with increased size, which affects the flow over the whole back end, causing a three-dimensional wake. These vortices are also responsible for maintaining attached flow over angled surface up to an angle of 30 • ; from 30 • and above, the flow is fully separated. There remains though a weak tendency of the flow to turn around the side edge of the model, a result of the relative separation positions of the flow over model top and that over the backlight side edges.
Due to some limitations on the wind tunnel and resources, Ahmed performed only force measurements on the bluff body during his experiments. In order to better understand the flow phenomena on an Ahmed Body, Lienhart and Becker [6] performed a study using Laser Doppler Anemometry (LDA), Hot-Wire Anemometry (HWA) and static pressure measurements in order to investigate the flow and turbulence structure around the Ahmed Body model for two slant angle conditions: 25 and 35 • . The main scope was to supply a detailed data set acquired under well-defined boundary conditions, similar to Ahmed first test, which considered a Reynolds number of 4.29 Million based on the length and static floor, to be used as reference data for numerical simulations.
Aiming to reproduce the real highway conditions of a vehicle, Strachan et al. [12], performed an Ahmed Body wind tunnel test with moving road conditions and both the aerodynamic forces and flow characteristics by time-averaged LDA were recorded. The flow conditions were slightly different from the ones used on Ahmed first test, by reducing the flow velocity to 25 m/s resulting in a Reynolds number of Re = 1.7 Million based on its length and the supports on the ground were replaced by a fixing system on the top of the tunnel, due to the rolling road simulation.
The Ahmed Body stands as one of the most used validation cases for CFD codes employed for automotive applications. Simulations employing a Reynolds Averaged Navier-Stokes (RANS) methodology are able to predict with good accuracy the drag coefficient, even for cases with complex flow topology, such as the slant angle of 25 • , with correlation factor of around 95% compared to experimental results, however the flow physics does not agree, usually under-estimating the flow features. Attempts considering more refined methodologies such as Detached Eddy Simulations (DES) and Large Eddies Simulation (LES) provide better correlation with experiments when comparing the flow structures but aerodynamic quantities values lose accuracy.
A trend that rose to improve the confidence level of CFD simulations was the high-order or high-fidelity methods, such as the spectral/hp element method [5]. The spectral/hp elemental method combines, according to Xu et al. [13], the advantages of the spectral element method, in terms of the properties of accuracy and rapid convergence, with those of the classical h-version finite element method, that allows complex geometries to be effectively captured. It also provides an attractive higherprecision approximation to solve partial differential equations.
One of the software that employs the spectral/hp element methodology is Nektar++ [9]. Nektar++ is a cross-platform spectral/hp element framework which aims to make high-order finite element methods accessible to the broader community. This is achieved by providing a structured hierarchy of C++ components, encapsulating the complexities of these methods, which can be readily applied to a range of application areas, as stated by Cantwell et al. [3]. It allows the use of high complex solution such as implicit LES (iLES) using a Spectral Vanish Viscosity (SVV) technique to stabilize the solution.
The latest achievements in the high-fidelity turbulence models around an Ahmed Body with slant angle of 25 • are summarized in the compilation work of Serre et al. [11], in which a comparative analysis of recent simulations, conducted in the framework of a French-German collaboration on LES of Complex Flows at Reynolds number of 768,000. It compares the results obtained with different eddyresolving modelling approaches, with two LES on body-fitted curvilinear grids: LES with Smagorinsky model and wall function (LES-NWM) and Wall-resolving LES with dynamic Smagorinsky model (LES-NWR), a stabilized spectral method known as iLES-SVV, similar to the one used in this present work, which is the base of the Nektar++ code and a DES-SST approach on an unstructured grid with element number ranging from 18.5 to 40 Million. Results of the flow field shows good agreement with results measured by Lienhart and Becker [6] by a gap on the drag coefficient values of 17% for the best case and 45% for the one using iLES-SVV.

Objectives
The main objective of this work is to evaluate the aerodynamic behaviour in terms of the drag and lift coefficients, considering an Ahmed Body with slant angle of 25 • using a spectral/hp elements method methodology as shown on Fig. 2. To achieve this, we first present a mesh study, evaluating two different size refinements referred as h-refinement for each of those, we employ three high-order surface mesh values to improve curvature representation. As the spectral/hp element method has also the possibility to improve the solution by increasing the polynomial order and consequently the number of degrees of freedom, we also evaluated three high polynomial orders for each mesh case, in a total of eighteen load cases.

Spectral/hp iLES-SVV
In this work, Nektar++ is used to run an implicit LES simulation using spectral/hp method. In this method, the domain is first divided into non-overlapping elements, offering geometric flexibility and allows for local refinement. Simulations were performed using the incompressible Navier-Stokes solver employing a velocity correction scheme, combined with a Continuous Galerkin (CG) projection. More details are presented by Cantwell et al. [3].
The mathematics behind Nektar++ basically considers the numerical solution of partial differential equations (PDEs) of on a domain , which may be geometrically complex, for some solution u. Practically, takes the form of a d-dimensional finite element mesh consisting of elements K i , embedded in a space of dimension dc, such that d ≤ dc ≤ 3, with = u i K i is an empty set or an interface between elements of dimension dbar < d. The PDE problem is solved then in the weak sense, considering that u|K i must be smooth with at least a 1st-order derivative. Therefore is required that u|K i is in the Sobolev space W 1,2 (K i ) equivalent to H 1 (K i ), according to Adams [1]. For a continuous discretisation, we impose C 0 continuity along element interfaces.
We assume the solution can be represented as u δ (x) = nû n n (x), a weighted sum of N trial functions n (x) defined on and the problem becomes that of finding the coefficientsû n . The approximation u δ does not directly give unique choices for the coefficientsû n . To achieve this, a restriction is placed on the residual so that its L2 inner product, with respect to the test functions n (x), is zero. For a Galerkin projection it is chosen that the test functions are the same as the trial functions, that is n = n . As outlined previously, to construct the global basis n it is first considered the contributions from each element in the domain. Each K i is mapped from a standard reference space K is between [−1, 1] by a parametric mapping χ e: K becomes K i given by x = χe(ξ ), where K is one of the supported region shapes, and ξ are d-dimensional coordinates representing positions in a reference element, distinguishing them from x which are d-dimensional coordinates in the Cartesian coordinate space.
The next step is to construct a local polynomial basis on each reference element with which to represent solutions. For 3D regions, a tensorial basis may be used, where the polynomial space is constructed as the tensor-product of one-dimensional bases on segments, quadrilaterals or hexahedral regions.
Spectral/hp element discretisation generally lead to approximations that have low dissipation and low dispersion per degree of freedom when compared to lowerorder methods. As stated by Xu et al. [13], in solving advection-diffusion equations and nonlinear partial differential equations such as advection-dominated flows, at marginal resolutions, oscillations appear that may render the computation unstable. Artificial viscosity has been used in may discretisation methods to suppress wiggles associated with high wavenumbers has been broadly and effectively used in simulations using the Fourier method. A related concept is the so-called SVV, which was originally proposed based on a second-order diffusion operator for spectral Fourier methods. SVV has been explicitly regarded as a turbulent model of implementing iLES under the assumption that the action of subgrid scales on the resolved scales is equivalent to strictly dissipative action stated by Sagaut [10], even though SVV is not explicitly designed as a subgrid-scale model. An example of a 1-D SVV kernel is: where P is the total number of modes employed and P cut is the cutoff polynomial order. SVV with the kernel function Df can be regarded as a low-pass filter. We see that the SVV dissipation added to the high mode numbers with respect to the spectral element discretisation does indeed yield dissipation at the global high wave number scales of the solution.
For this work, we employed a novel CG-SVV scheme with DGKernel, proposed by Moura et al. [8] where he dissipation curves of CG of order p are match to those of DG with order p − 2, eliminating non-smooth dissipation characteristics arising from CG dissipation when considering high Reynolds number.

Simulation Methodology
We first define the coordinate system as X the streamwise direction, Y the vertical direction and Z the spamwise direction. The Ahmed Body length of 1.044 m is defined as 1 AL. The virtual wind tunnel dimensions are 2.74 × 1.66 m for the test section and total length of the domain of 4 AL, similar to Strachan et al. [12] study. The Ahmed Body model back in placed on X = 0, inlet position at X = −2 AL and outlet position at X = 2 AL. A schematic setup is shown on Fig. 3.
In terms of boundary conditions, velocity was normalized to 1 in order to match the Reynolds number previously stated and set as the inlet boundary condition. The outlet was set as pressure high-order outlet condition and the floor was also set with the same velocity of the free stream in order to reproduce the moving floor effect. The top and outer side wall and the Ahmed Body wall are set as no slip condition and a symmetry condition. Total simulated time is 7 convective lengths AL, which means that the flow is able to cross the whole domain.
This study evaluates two mesh configuration considering different h-refinements and referred as Original and Refined meshes and for each of those, three highorder surface mesh settings: 4th, 5th, and 6th order, generating six different meshes. All mesh files were generated by NekMesh, which is Nektar++ high-order mesh generator. In both Original and Refined meshes, cases two refinement zones were generated, where the first one, defined as the Ahmed Body refinement, ranges from 0.3 AL before the beginning of the geometry and 0.3 AL after the end of the body, in a total length of 1.6 AL. The second refinement, defined as the Wake Refinement region, intercepting the first refinement in 0.3 AL before the end of the body, to 1.3 AL after the end of the body, in order to fully capture the flow phenomena in the separation region, with same total length of 1.6 AL, as illustrated on Fig. 4.
The Original Mesh has total number of elements for half model around 95,000. For the Refined mesh, the boundary layer setup was the same and the dimensions   Most of the commercial CFD code employ low order methods and the highest order polynomial interpolation for the solutions usually seem is 3rd. The mesh plays the major role for complex simulations such as LES, leading to elevate number of elements to reach a reliable result. To make use of the flexibility of the spectral/hp element methods, we proposed solutions considering polynomials with order higher than 3rd within the previous mesh refinement studies as the higher order polynomials increase the degrees of freedom and resolution of the mesh. For the Nektar++ implicit LES simulations using the Incompressible Navier-Stokes solver evaluated three different polynomial expansions, 4th, 5th and 6th orders, referred here as P4, P5 and P6. In summary, 18 load cases were evaluated using HPC with 432 CPUs for each case.

Drag Coefficient Comparison Results
The drag coefficient for the 18 cases evaluated, considering 9 from the Original mesh with 95,000 elements considering fourth, fifth and sixth polynomial order and the Refined mesh case with 310,000 elements also considering fourth, fifth and sixth polynomial order expansion, with maximum RMS and compared with experimental results are shown on Fig. 6.
From Fig. 6 it is possible to observe that for the drag coefficient, P4 polynomial expansion considering both mesh cases presented mean drag results around 35% higher than the experimental results. For the P5 cases, considering again both Original and Refined mesh cases, the error was reduced to 5% however results change the trend from over-predicted to under-predicted when the mesh is refined further. The cases considering P6 polynomial expansion presented the same trend for both mesh cases, highlighting its consistency although the mean error when compared to experiments increases to 16%.

Lift Coefficient Comparison Results
Similar to the drag coefficient graph, in Fig. 7 the lift coefficient for the all evaluated cases is shown, considering Original and Refined meshes and fourth, fifth and sixth polynomial order expansion. Maximum RMS is also plotted for all cases and compared with experimental results from Strachan et al. [12].  Analyzing Fig. 7, we observe that the h-refinement from Original mesh to Refined mesh lead to results closer to experimental values when adopting P4 as the polynomial expansion basis. For both P5 and P6 polynomial expansions, lift coefficient results present good agreement with experimental data, with maximum mean error of 5%.

Flow Structure Comparison
In terms of the polynomial order expansions for the solution, combined with the 6th order surface mesh, the results present focus on the Refined mesh case, once they improved the correlation for the P4 polynomial expansion within experimental results and kept similar trend for P5 and consistent results for P6 in terms of drag and lift coefficient prediction. An initial comparison in terms of flow structures is present by the Q-Criterion of 350 coloured by pressure, comparing the Refined mesh case, considering P4, P5 and P6 polynomial expansions in Fig. 8.  Fig. 8 it is possible to visualize that P4 is unable to define the vortex on the side of the slant, explaining also the difference in terms of both drag and lift coefficients, compared to experimental results. Results for P5 show the side vortex clearly defined and P6 is also able to capture the lower vortex, detailed on the lower image, which is not present in the studies considering the Ahmed Body, but they are important to understand the behaviour with the moving floor. Figure 9 shows a contour of Lambda 2 to illustrate the lower vortex detail on the plane x/L = 0, on the back of the Ahmed Body.
We next focus only on the Refined mesh with 6th order surface mesh for P5 and P6, once they were able to predict both lower and top vortices. Due to nature of the wind tunnel with moving ground used by Strachan et al. [12], the model had to be fixed on the top by a steam, which can be removed in the drag coefficient calculations, however it might change the flow topology over the slant, as stated by the authors themselves.
Comparing the plane x/L = 0.076 with the measurements of the flow velocity on x direction U normalized by the free stream velocity of Lienhart and Becker [6] with static floor without the steam on the upper portion with results of Strachan et al. [12], it is possible to notice intensity changes in the U normalized velocity and this is attributed by the last due to the upper support.
As the simulations do not included the upper support, but do include the moving ground, the expected results are the top portion to be similar to Lienhart and Becker [6] measurements and the lower part, correlated with measurements of Strachan et al. [12], which both P5 and P6, proved to have good agreement in terms of normalized U velocity, shown in Fig. 10.
Similar comparison is presented on Fig. 11 for the vortex intensity on the slant, on the plane x/L = 0, on the back of the Ahmed Body for vertical velocity V normalized by the free stream velocity. In this, the simulations close correlate to Lienhart and Becker [6] study, due to the absence of the steam support but for this case, the higher polynomial order expansion P6 is able to capture more scales than the P5 for the core of the main vortex, highlighting the gain of resolution of the high-order simulations.  [6] with static floor without the steam (left), results of Strachan et al. [12] with moving floor and steam support (middle left), Refined mesh and 6th order surface P5 (middle right) and Refined mesh and 6th order surface P6 (right)

Conclusions
Within the advances in CFD codes, confidence level and computational power, aerodynamic simulations are applied in almost every automotive company. The reason is very simple: reduced development cost and time, which is an enormous advantage in a competitive market. High-fidelity simulations are becoming a reality for complex industrial cases in order to improve resolution and results in a reliable response time, such as presented for the Ahmed Body on this work.
On the meshing definition study, the surface mesh order seems not to influence the results in terms of aerodynamic quantities, presenting similar trend for same polynomial order, as the Ahmed Body geometry has curved surfaces only on the front portion.
Still on the mesh definition, as the h-refinement increases from Original to Refined mesh, the drag coefficient values for P4 and P6 remains unchanged and P5 values switched from positive to negative. We conclude that consistency is shown for P4 and P6 cases but P6 presented the most reliable results, with a maximum deviation of 16%. For the lift coefficient, results for P4 improved as the h-refinement increased and kept similar values for both P5 and P6, where the best agreement was found for the case considering Refined mesh with 6th order surface mesh and P6 as the polynomial expansion.
Flow structure results focus only the Refined 6th order surface mesh, where the main expected features were captured by P5 and P6 cases. It was confirmed by those two simulation cases that the lower portion has similar behaviour of the moving ground test conducted by Strachan et al. [12], however the top portion close correlate to Lienhart and Becker [6] experiments, as the simulation cases allow the body to be fixed without the upper support used in the experiment. This fact might also explain the difference from the simulation results with the literature experiments, as the simulation allows idealized configurations.
For all simulation cases, half of the body is being simulated and a symmetry plane is set on the middle portion. From Fig. 12, which shows the normalized U velocity on a line of coordinate y/L=0.15 on the plane x/L = 0, we observe that simulation has good agreement with experimental results, with a small distortion as it gets closer to the symmetry plane.
Open Access This chapter is licensed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence and indicate if changes were made. The images or other third party material in this chapter are included in the chapter's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the chapter's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.