12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022

PHYSICS BASED TAILORED NUMERICAL GRID FOR STATISTICALTURBULENCE MODELS BASED ON DISSIPATION ELEMENTS
Fettah Aldudak, Holger Foysi
Department of Mechanical Engineering, University of SiegenPaul-Bonatz-Straße 9-11, 57076 Siegen, Germany
fettah.aldudak@uni-siegen.de, holger.foysi@uni-siegen.de
ABSTRACT

cending and descending scalar gradients will always reach aminimum and a maximum point. A DE is defined by all points
along the trajectories belonging to the same pair of extremal
points. Hence, DEs directly originate from the turbulent flow
and, at the same time, are irregularly shaped and space-filling.
This means that a turbulent scalar field can be completely
decomposed into such elements with no empty space in between. Velocity components, turbulent kinetic energy k andits dissipation rate ε can be chosen for such a decomposition.
The DE concept has been introduced first by Peters & Wang
(2006) who analysed small-scale statistics from passive scalar
fields of homogeneous shear turbulence obtained from DNS.
For this flow case, Wang & Peters (2006) report that the mean
DE length ℓDE is of the order of the Taylor scale defined as
λ = (10νk/ε)1/2 where ν is the kinematic viscosity, k the turbulent kinetic energy and ε the dissipation rate. As the DE
length scales cannot be derived accurately when using RANS,
for example, a relationship like this would allow us to base thelength scale determination on a quantity easily accessible evenin statistical simulations.
The DE methodology has been further applied within
DNS of canonical turbulent channel flow by companion worksof Aldudak & Oberlack (2009), Aldudak & Oberlack (2011),
Aldudak & Oberlack (2012) and Aldudak (2012). The strong
influence of rigid walls on the flow geometry and structures
has been examined showing a distinct dependency of the DE
length on the wall-normal direction. It was shown that ℓDE
increases linearly in the wall-normal direction y, starting fromthe logarithmic region, until a constant length is reached in
the channel center, due to decreasing shear. Interestingly, thislinear behaviour occurs mainly in the intermediate region between the inner and outer layers, where the turbulent production of k and its dissipation ε are approximately equal at sufficiently high Reynolds numbers (Hinze, 1975). This regionis roughly located between y+ = uτ y/ν = 30 and y/h = 0.7(Aldudak, 2012). This is in good agreement with the findings in Lee & Moser (2015) where the balance of production and dissipation of k in a turbulent channel flow DNS at
Reτ = uτ h/ν = 5200 is plotted against the wall-normal
direcption. Here, h is the channel half width, uτ = τ/ρ is thefriction velocity and τ and ρ are the mean wall friction anddensity, respectively. The region of relative balance coincides
well with the above-mentioned range of linear increase of themean DE length scale. Furthermore, it is found that the DE
length scale decreases with increasing Reynolds number Re,
i.e. ℓDE = f (y, Re), similar to the Taylor length λ . Given thespace-filling nature of the DEs, this means that the numberof generated structures in a wall-bounded turbulent flow decreases with wall distance and increases with Reynolds number.
In Aldudak (2012), we investigated the relation between
DE and classical length scales, especially the intermediate

Since direct numerical simulations (DNS) remain out of reachfor most realistic flow conditions, statistical methods as employed in Reynolds Averaged Navier Stokes (RANS), Detached Eddy Simulation (DES) or Large Eddy Simulation(LES) represent an alternative. Regardless of the choice, quality and resolution of the corresponding numerical grid are crucial for the accuracy of a simulation. In the case of DNS, there
exists a very well defined requirement in the form of the Kolmogorov length scale η to guide the numerical grid resolution,for example. Unfortunately, for statistical turbulence models,finding the right resolution criterion for automatic grid refinement is more difficult. In addition, real flow conditions require
locally adapted, inhomogeneous grids to account for complex
flow regimes. In this work, we propose the application of theso-called Dissipation Element (DE) scale ℓDE and its relation
to the Taylor length scale λ as a rough guide for the grid resolution for the use in automatic grid adaptation methods. Both,canonical and realistic complex flows are used to propose agenerally applicable condition for the choice of the grid resolution. The present method automatically ensures an adapted
resolution throughout the domain with local refinements, minimizing necessary user input to a single quantity.

INTRODUCTION
Generating a high-quality numerical grid is an integral partof the simulation workflow which usually requires a significant amount of time. To ensure grid independent numericalsolutions, several cycles of grid generation and subsequent
flow simulation are performed iteratively until no substantialchanges in the results are observed. Therefore, shortening thiseffort is accompanied by a reduction of the temporal and computational costs.
The presented method aims at obtaining the targeted grid
faster while additionally allowing the appropriate resolution
locally. The latter is essential for the efficiency since refinement is performed only where necessary. In a similar manner to DNS resolution down to the Kolmogorov length η, we
intend to find a physically-based length scale for turbulencemodels to adapt the grid resolution based on physics. Ideally,a universally valid magnitude of this length scale should directly emanate from the flow physics without the need of further user input. Thus, adapting the resolution to this global
value will lead to the desired unique mesh for the flow underconsideration which is ready to use for the final simulation.
We propose the Dissipation Element (DE) to be such aunique, physically-based length scale. These intermediate spatial structures associated with turbulent scalar fields can be
identified by connecting corresponding pairs of zero-gradientextremal points. Finite-length gradient trajectories startingfrom every point in the scalar field in the directions of as-1

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022
A surface mounted cube (Martinuzzi & Tropea (1993))
was simulated using RANS, where no-slip boundary condition was applied to the horizontal top, bottom walls and the
cube, respectively. The Reynolds number based on the cube
height H was chosen to be ReH = 40000. The spatial dimensions in the streamwise, spanwise and wall-normal directions
were 19.5H × 9H × 2H. Further information can be found in
Martinuzzi & Tropea (1993) where the flow was investigatedexperimentally. Figure 2(i) depicts the coarse initial grid using uniform cubic cells starting the simulation. Figure 2(ii)instead shows the final grid where all cells were refined to fulfill the resolution criterion in (1). Near the walls and aroundthe cube, as well as in regions with large gradient, the refinement level is highest, transitioning to less refined regions using consecutive refinement levels, as expected. Figure 2(iii)
shows results of the velocity obtained from the unsteady k-ωSST RANS (Langtry & Menter, 2009) and Spalart-AllmarasIDDES DES (Gritskevich et al., 2012) simulations performed
with OpenFOAM and our grid refinement procedure, demonstrating good agreement with experimental data from Martinuzzi & Tropea (1993).
In a similar manner, flow past an Ahmed body (Ahmed et al.
(1984)) configuration with a 25◦ slanted back face was calculated with k-ω-SST RANS which, in spite of its simple geometry, exhibits complex flow aspects associated with car aerodynamics such as a counter-rotating vortex pair and a strong separation and recirculation region. The Reynolds number basedon the length of the body was Re = 2.8 × 106 . Figure 3(i)shows two recirculation bubbles, as expected, together withthe case-specific tailored mesh and the corresponding velocity
profiles compared with experimental data. The tailored gridshows refinement regions of different intensity according to
the resolution criterion in Equation 1 with finer cells near walls
as well as separation zones. The typical longitudinal counterrotating vortices (Lienhart et al. (2002)) in the body wake areaccounted for by the grid refinement process. In the last figure, mean velocity profiles obtained from RANS simulations
are plotted together with measurements from Lienhart et al.
(2002) at selected positions along the symmetry plane on theslant and in the wake of the Ahmed body. The profiles are ingood agreements with the experiments, especially in the wakeregions with recirculation, whereas near the trailing edge ofthe slant small deviations are visible.

Taylor length. It was found that the DE length is proportional
to the Taylor scale in turbulent channel flows almost overthe entire channel height (ℓDE ∼ λ ). To find an expectedvalue for the proportionality factor, we plot the ratio ℓDE /λas a function of wall-normal distance for the friction based
Reynolds number Reτ = 720. The three plots represent theratio between DE and Taylor scales for individual velocity
components, respectively. As seen in figure 1, ℓDE /λ doesn’tvary widely after the start of the logarithmic region at abouty/h ≈ 0.04, with the proportionality factor being between
2 and 3 for a large wall-normal extent. An increase can beobserved close to the wall towards the start of the logarithmic
and subsequently buffer layer, where the production peaks,
too. Although each length scale itself exhibits a clear ydependency, their quotient appears to be widely independentof the wall-normal direction throughout a large part of thechannel height. This observation is a direct indication that thedissipation element length scale can be determined to good
approximation by making use of the Taylor microscale, whichis readily available even in Reynolds averaged simulations, forexample. Since the grids were often including points within
the buffer layer, too, the criterion for grid refinement waschosen to require cells sizes following the inequality∆ < ℓDE ≈ 3 · λ ,(1)

to reflect the increase of ℓDE /λ close to the wall. This inequality is not exact but turned out to give nice results in various
numerical experiments.
To verify the ability of the proposed method several
canonical and more realistic test cases were simulated using
RANS and DES. So far, the test cases included plane turbulent channel flow (not shown here), a surface-mounted cubein a channel (Martinuzzi & Tropea (1993)), an Ahmed body
generic car configuration (Ahmed et al. (1984)) and a more demanding case of a passenger car (DrivAer, Heft et al. (2012a)).
The latter cases exhibit complex flow structures with large separation and recirculation regions as well as sudden changes ingradients. At the same time, they are also well-known benchmarks, thus suitable for validation purposes.

RESULTSA solver in the framework of OpenFOAM (Weller et al.
(1998)) was extended to identify and automatically refine numerical cells during the simulation process, whose size violated inequality (1).

Since the two aforementioned geometries are rather simple, the more demanding case of a passenger car has been used
to further validate the proposed tailored mesh approach. With
this test case we intend to demonstrate the efficiency of thisnew method further for relevant, complex flow problems. As aconfiguration the realistic DrivAer car model is chosen which
was introduced by Heft et al. (2012b). Its modular geometry
with three different rear end configurations, i.e. estate back,notchback, and fastback, allows for aerodynamic research on
realistic passenger cars. Details of varying degree, such as mirrors, wheel rims and a detailed underbody can be added. Different combinations of DrivAer models have been investigatedboth numerically and experimentally making it eligible for validation studies. For the present work, the fastback variant withsmooth underbody, side mirrors and detailed wheels is utilized
and results are compared with available data. As before, an
initial numerical grid is generated first which is adapted tothe shape of the car and consists of hexaeder-dominant coarse
cells. A precursor simulation has been conducted to provide
appropriate turbulent inflow profiles. Subsequently, a first simulation has been performed, yielding the cell to Taylor length
ratio, ∆/λ , as the key parameter for local refinement. Based on(i) Initially, a very coarse mesh with uniform grid cell size wascreated to serve as a starting numerical grid. The cells are cubic so that the length of the cell, ∆, is calculated using the cubicroot of the individual cell volume.
(ii) After the simulation reached an intermediate steady state,
the algorithm is refining all those cells further, which exceed
the resolution criterion. Cells smaller than half the resolutioncriterion are coarsened.
(iii) After every refinement step, the old solution is interpolated linearly onto the new mesh by the solver and the computation is resumed until a new intermediate steady state isreached.
Step (ii) and (iii) are repeated, until all cells fulfill the resolution criterion and steady state is reached.
The algorithm was capable of successfully producing readyfor-use grids with well-distributed local refinement regions using inherent flow physics only.

2

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022
tal data until the trailing edge of the rear end. While the LES
results match the measurements more accurately, minor deviations in RANS are visible. Similarly, figure 7(b) illustrates thepressure distribution at the bottom of the symmetry plane. A
pronounced difference between the measurements and simulation data in general is observed initially, where the measured
C p decreases sharply. Similarly to the top of the car, the RANS
results exhibit an overall good agreement with LES and experimental data, validating the automatic grid refinement procedure.

this quantity, all cells exceeding the criterion in Equation 1 are
identified and refined accordingly. Thereby, the desired final
mesh with appropriate resolution is obtained and can be used
for further simulations. The Reynolds number based on the
car length is set to 4.87 × 106 which was investigated experimentally by Heft et al. (2012b), Strangfeld et al. (2013)) andnumerically by Rüttgers et al. (2019), Ekman et al. (2020), forexample. The final mesh for a symmetric half-body contains
53 Million cells with the first wall boundary layer located in
the lower logarithmic layer, hence, giving rise to the use ofstandard wall functions. The size of the mesh is comparable
to that used in a Large-Eddy-Simulation with 55 Million cells,
discussed in Rüttgers et al. (2019). The mean inflow velocity is 16 m/s and a no-slip boundary condition is applied on
rigid walls. The two-equation k − ω−SST turbulence model
has been employed for a steady-state RANS simulation.
Figure 4 shows iso-contours of two regions with the highest refinement levels which are identified by the tailored meshapproach. The purple region around the side mirrors and thenear-wake marks the finest grid resolution and is located in therecirculation region which normally forms behind bluff bodies.
Mirrors, wheels and the rear end are zones where elongatedvortical structures are expected, see Rüttgers et al. (2019). The
larger gray iso-surface starts from the wheels and the side mirrors and encloses large areas around the car. Note, that typical counter-rotating longitudinal vortex regions in the wake
are detected for proper resolution, illustrating the physicallybased nature of the proposed method. As before, individualrefinement regions are created with smooth transitions as seenin Figure 5 and the zoomed region behind the car in the same
figure.
In image 6, streamlines based on the averaged velocity
on the symmetry plane are depicted together with the wakeisoline (black) where the streamwise velocity is zero, highlighting the common recirculation zone in the near-wake(Rüttgers et al. (2019), Ekman et al. (2020), Strangfeld et al.
(2013)). Moreover, streamlines coming from the upper part
form a small separation zone towards the end of the rearwindow. Additionally, the pressure coefficient defined asp − p∞
Cp =0.5ρU∞2

CONCLUSIONS
A novel tailored grid approach has been introduced,
where Dissipation Element (DE) length scales ℓDE , obtainedby DNS results, serve as grid resolution criterion. Previous
work has been discussed to underline the relation between DE
and the well-known Taylor length scale (λ ) which is usuallyavailable to statistical simulation methods. Based on resultsfrom turbulent channel flow simulations, where both lengthscales can be accurately determined, it was found that the relation of DE and Taylor length scale ℓDE /λ is varying only toa minor degree over most of the channel width, starting withthe logarithmic region. This allows determining ℓDE via the
Taylor scale to automatically refine the grid whenever the gridcells surpass a predetermined threshold based on ℓDE .
Benchmarks and test cases such as a wall-mounted cube,
Ahmed body and a more realistic DrivAer passenger car havebeen analyzed by applying the proposed meshing method. In
all cases it could be observed that, by applying the resolution
criterion, a well-formed numerical grid with smooth transition
emerged with no further involvement required externally. The
resulting numerical grids are shown to follow physical characteristics of the underlying turbulent flow closely and are capable of reproducing the validation results with appropriate accuracy. Contrary to classical meshing techniques, expensive gridindependence studies can be reduced to a minimum. Since
scale-resolving simulations tend to DNS with increasing resolution, this proposed resolution criterion could serve as a goodcompromise between direct resolution and modelling. Nevertheless, further analysis is needed to investigate the tailoredgrid approach in conjunction with LES and DES.

is plotted on the car surface where U∞ , p∞ , and ρ are inflow
velocity, freestream pressure and and air density, respectively.
In order to look at the pressure distribution in more detail andcompare to the literature, C p is shown in figure 7 for various
streamwise positions along the car determined at the symmetry
plane. Experimental data from Heft et al. (2012b) and numerical results of a LES by Rüttgers et al. (2019) are comparedhere to the present RANS results obtained with the tailoredgrid method. Figure 7(a) shows the distribution of the pressure coefficient C p along the symmetry line of the DrivAer car
at the top surface. A good overall agreement of our RANS
results with both experiments and LES data can be seen. At
the stagnation point at the vehicle front a high pressure isobserved, decreasing rapidly first due to the acceleration of
the flow. Along the adjacent engine bonnet, the pressure increase is replicated appropriately. Additionally, the distribution around the cowl top, showing a small kink in the C p distribution around x/L = 0.3, is reproduced nicely with good accuracy, where separation and reattachment occur. Arriving at
the car roof, both LES and RANS deviate by the sudden pressure change and underpredict the pressure. To some extent,
this deviation can be caused by the support holding the car in
the experimental installation (Heft et al. (2012a), Ashton et al.
(2016)). The pressure distribution is close to the experimen-

REFERENCES
Ahmed, S.R., Ramm, G. & Faltin, G. 1984 Some salient features of the time-averaged ground vehicle wake. In SAE
International Congress and Exposition. SAE International.
Aldudak, F. 2012 Geometrical structure of small scales
and wall-bounded turbulence. PhD thesis, TU Darmstadt,
Darmstadt, Germany.
Aldudak, F. & Oberlack, M. 2009 Dissipation element analysis
of scalar fields in wall-bounded turbulent flow. Proceedingsof EUROMECH Colloquium 512, pp. 9–11.
Aldudak, Fettah & Oberlack, Martin 2011 Length scale analysis in wall-bounded turbulent flow by means of dissipation
elements. J. of Physics: Conf. Series 318 (2), 022001.
Aldudak, F. & Oberlack, M. 2012 Dissipation element analysis
in turbulent channel flow. J. Fluid Mech. 694, 332–351.
Ashton, N., West, A., Lardeau, S. & Revell, A. 2016 Assessment of rans and des methods for realistic automotive models. Computers Fluids 128, 1–15.
Ekman, P., Wieser, D., Virdung, T. & Karlsson, M. 2020 Assessment of hybrid rans-les methods for accurate automotive aerodynamic simulations. Journal of Wind Engineering
and Industrial Aerodynamics 206, 104301.

3

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022
Gritskevich, M.S., Garbaruk, A.V. & Schütze, J. 2012 Development of ddes and iddes formulations for the k − ωshear stress transport model. Flow Turbulence Combust. 88,
431–449.
Heft, Angelina I., Indinger, Thomas & Adams, Nikolaus A.
2012a Experimental and numerical investigation of the drivaer model. Fluids Engineering Division Summer Meeting,vol. Volume 1: Symposia, Parts A and B, pp. 41–51.
Heft, Angelina I., Indinger, Thomas & Adams, Nikolaus A.
2012b Introduction of a new realistic generic car model for
aerodynamic investigations. In SAE 2012 World Congress
Exhibition. SAE International.
Hinze, J. O. 1975 Turbulence. McGraw-Hill Book Company.
Langtry, Robin B. & Menter, Florian R. 2009 Correlationbased transition modeling for unstructured parallelized
computational fluid dynamics codes. AIAA Journal 47 (12),2894–2906.
Lee, Myoungkyu & Moser, Robert D. 2015 Direct numerical simulation of turbulent channel flow up to reτ =5200.
Journal of Fluid Mechanics 774, 395–415.
Lienhart, Hermann, Stoots, C & Becker, S. 2002 Flow andturbulence structures in the wake of a simplified car model(ahmed model). Notes on Numerical Fluid Mechanics 77.
Martinuzzi, R. & Tropea, C. 1993 The Flow Around SurfaceMounted, Prismatic Obstacles Placed in a Fully Developed
Channel Flow. J. of Fluids Engineering 115 (1), 85–92.
Peters, N. & Wang, L. 2006 Dissipation element analysis ofscalar fields in turbulence. C. R. Mech. 334, 493–506.
Rüttgers, Mario, Park, Junshin & You, Donghyun 2019 Largeeddy simulation of turbulent flow over the drivaer fastback
vehicle model. Journal of Wind Engineering and Industrial
Aerodynamics 186, 123–138.
Strangfeld, Christoph, Wieser, Dirk, Schmidt, Hanns-Joachim,
Woszidlo, Rene, Nayeri, Christian & Paschereit, Christian
2013 Experimental study of baseline flow characteristics
for the realistic car model drivaer. In SAE 2013 World
Congress Exhibition. SAE International.
Wang, L. & Peters, N. 2006 The length scale distributionfunction of the distance between extremal points in passive
scalar turbulence. J. Fluid Mech. 554, 457–475.
Weller, H. G., Tabor, G., Jasak, H. & Fureby, C. 1998 A tensorial approach to computational continuum mechanics using
object-oriented techniques.

(i)

(ii)

z/H

(iii)

x/H

Figure 2: Wall mounted cube: (i) Starting grid. (ii) Tailoredgrid after refinement. (iii) Velocity profiles at selected positions, comparison to experimental data of Martinuzzi & Tropea (1993).
(i)ℓDE /λ

(ii)

z/H

(iii)

y/h

Figure 1: Proportionality relation between DE and classicalTaylor length scales for Reτ = 720 over wall-normal distancey/h obtained from DNS results (Aldudak, 2012). Gray area
starts at the log-layer and marks the broad region of validity
for the refinement condition.

x/H

Figure 3: (i) Streamlines at the center plane. (ii) Tailored grid
after automatic refinement. (iii) Velocity profiles at selected
positions in comparison to data from Lienhart et al. (2002).

4

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022

Figure 4: Illustration of the highest refinement regions. The light grey colors represent the second finest grid region, the purple region
the finest grid cells. Those are located behind the side mirrors and tail of the car. Contours of C p are depicted on the car surface.

Figure 5: Tailored grid of the fastback DrivAer configuration (Heft et al. (2012a)) in the symmetry plane, showing different levels ofrefinement. Top: Tailored grid, complete car. Bottom: Zoom into the tail section of the car.

5

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan, July 19–22, 2022

Figure 6: Surface pressure coefficient distribution, mean velocity streamlines along the symmetry plane and wake isoline Ux = 0(black).

(a) Pressure coefficient at the top surface.

(b) Pressure coefficient at the bottom surface.
Figure 7: Pressure coefficient distributions along the symmetry plane of the DrivAer car (indicated in grey), for the top (a) and bottom(b) surface, compared to experiments (Heft et al. (2012a)) and LES simulations (Rüttgers et al. (2019)). L indicates the length of thecar.

6