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

A SUB-GRID ACTIVITY SENSOR APPLIED TO MIXED MODELING IN
LARGE EDDY SIMULATION

Josef Hasslberger
Institute of Applied Mathematics and Scientific Computing
Department of Aerospace Engineering, University of the Bundeswehr Munich
Werner-Heisenberg-Weg 39, 85577 Neubiberg, Germany
josef.hasslberger@unibw.de

SUB-GRID SCALE MODELINGIn contrast to functional eddy viscosity models, structuralmodels aim to reproduce the structure of the SGS stress tensor itself and, among those, the Bardina/Liu model (Bardina,
1983; Liu et al., 1994)

Mixed models, i.e. a combination of functional and structural sub-grid scale (SGS) models, are among the most successful approaches for Large Eddy Simulation (LES). Compared to purely structural models, mixed models are especiallysuperior in terms of stability. Chapelier et al. (2018) have recently demonstrated the potential of a sub-grid activity sensor
to improve the performance of functional eddy viscosity models in regions with transitional flow features. Further, it has
been shown by Hasslberger et al. (2021) how a sub-grid activity sensor can additionally be used to rectify the incorrect nearwall scaling of eddy viscosity base models like the standard
Smagorinsky model without explicit wall damping. Accordingly, the idea here is to exploit the advantages of a sub-gridactivity sensor in the context of mixed SGS modeling.
The coherent structure function, as required in the following analysis, is a useful quantity to characterize the structureof turbulent flows. It is defined as FCS = Q/E, i.e. the second
invariant of the grid-scale velocity gradient tensor
Q = (Ωi j Ωi j − Si j Si j )/2

d bb
τiSS
j = ui u j − ui u j

is based on the scale similarity (SS) principle, where (b·) represents a suitably defined explicit test filter. In this work, theexplicit test filter for any field quantity φi, j,k at the discrete
location given by the index triple (i, j, k) is implemented according to Anderson & Domaradzki (2012):+1

φbi, j,k = ∑q

2Skl Skl Si j = −2νt Si j

+1

(5)

This three-dimensional filter is the product of the convolution of three one-dimensional filters with coefficients(b−1 , b0 , b+1 ) = (c f il , 1 − 2c f il , c f il ) where c f il = 1/12.
The following blending scheme is based on the below observation that structural models of the scale similarity (SS)type perform the best for moderate under-resolution of theflow, whereas the concept of eddy viscosity (EV) becomes increasingly valid the lower the relative resolution is. Hence, theidea is to retain the superior alignment properties of SS models
in those regions where an extrapolation based on the smallest
resolved scales (which are also strongly affected by numericalerrors) is still properly working, and to blend in a more robust
EV model where physically plausible. Using a local sub-gridactivity sensor Θ, the blending scheme reads

(1)

(2)

where Si j = (∂ ui /∂ x j + ∂ u j /∂ xi )/2 is the grid-scale strain
tensor and Ωi j = (∂ ui /∂ x j − ∂ u j /∂ xi )/2 is the grid-scale rotation tensor. Consequently, FCS exhibits a definite lower andupper limit, i.e. −1 ≤ FCS ≤ +1. The values −1 and +1 correspond to pure elongation/strain and pure rotation, respectively. It is important to note that the constituents of Q and
E are related to fundamental quantities to describe turbulent
flows, namely the dissipation of kinetic energy into heat perunit mass, 2νSi j Si j , and enstrophy, Ωi j Ωi j .
The best-known functional SGS model is the standard
Smagorinsky model (Smagorinsky, 1963) which calculates the1
deviatoric part τidev
j = τi j − 3 τkk δi j of the SGS stress tensor forincompressible flows τi j = ui u j − ui u j as22
τiEV
j = −2 CSmago ∆+1∑ ∑ bl · bm · bn · φi+l, j+m,k+nl=−1 m=−1 n=−1

being normalized by its magnitude
E = (Ωi j Ωi j + Si j Si j )/2(4)SS= τiEVτimixed,sensor
j Θ + τi j (1 − Θ)
j

(6)

SS
where τiEV
j and τi j represent any kind of eddy viscosity andscale similarity type model, respectively. The lower the relative resolution, the higher the sub-grid activity Θ and thehigher (lower) the EV (SS) contribution is. This is also consistent from a physical point of view, because “randomly” fluctuating incoherent turbulence acts like diffusive motion – agreeing with the way eddy viscosity is reflected in the diffusive
term of the filtered Navier-Stokes equations.

(3)

where ∆ is the grid size (i.e. the implicit filter width), νt is theeddy viscosity (EV) and CSmago = 0.17 is the theoretical valueof the Smagorinsky constant.

1

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19–22, 2022
flows. Also, the SGS energy transfer ε = τi j Si j can be decomposed into forward scatter ε − = 0.5(ε −|ε|) and backward
scatter ε + = 0.5(ε + |ε|). For consistency, only the deviatoricpart of SGS stress is considered for the EV model, although1
the discrepancy between τidev
j = τi j − 3 τkk δi j and τi j is small
here.
Figure 3 shows that for moderate filter width (∆/∆DNS =9), the Bardina/Liu model clearly outperforms the Smagorinsky model, because the latter is overly dissipative. For large
filter width (∆/∆DNS = 33), both types of models are in good
agreement with the reference data. Under these conditions, the
SGS energy transfer, on average, is almost linearly increasingfrom rotation- to strain-dominated regions and this can be represented also by the Smagorinsky model. Although not explicitly shown here, tendencies with respect to Reynolds number variation are expected to be similar than for filter widthvariation. It cannot be seen from this a-priori analysis but it
is known that structural models tend to become unstable under high filter width/Reynolds number conditions which motivates the mixed model. It is worth noting that the SS model isable to represent backward scatter in contrast to the EV model.
Backward scatter is not unphysical but is discussed as a potential source of instability for structural models in the literature(Kobayashi, 2018; Klein et al., 2020).

The bounded sensor function Θ is constructed as

1




σ < σeq1 
σeq − 2σ + 1
1
+
sin
π
Θ=
if σ ∈ [σeq , 1] (7)22(1 − σeq )

σ >1

0
such that a smooth transition between well-resolved (σ > 1)
and insufficiently-resolved (σ < σeq ) regions is obtained. Although Eq. 7 is the same as in (Chapelier et al., 2018), σ iscalculated in a different manner. Rather than using enstrophyb
only, σ = E/E
accounts for both dissipation and enstrophy.
Comparison of the implicitly grid-filtered value E and the explicitly test-filtered value Eb allows to estimate the local subgrid activity. At the same time, the intensity-preserving (E ≈const.; Eb ≈ const.) natural exchange between dissipation andenstrophy in turbulent flows remains undetected by the sensor.
According to the definition of E, Eq. 2, fluctuations in strainand rotation-dominated flow regions can be equally well reflected by the sensor. This is demonstrated in Fig. 1 which
depicts conditionally averaged values of the sub-grid activity
Θ in homogeneous isotropic turbulence (details on this a-priori
analysis are provided subsequently). The original sensor formulation based on enstrophy only (left panel of Fig. 1) is obviously unable to detect sub-grid activity in strain-dominatedflow regions. Independent of the implicit filter width in LES,
strain-dominated flow regions are even more probable thanrotation-dominated flow regions in turbulent flows as shown by
the probability density function (PDF) of the coherent structure function FCS in Fig. 2. Agreeing with expectations, Fig.
1 demonstrates increasing levels of sub-grid activity Θ for increasing under-resolution of the flow as specified by ∆/∆DNS .
The calculation of the equilibrium value σeq is identical to Chapelier et al. (2018) since dissipation and enstrophy are obeying the same spectral scaling. Both spectra areproportional to κ 2 e(κ) with wavenumber κ and energy spectrum e(κ). Their peak is located at high wavenumbers (small
scales), hence also the naming as a sub-grid activity sensor.
On average, Eb < E in the inertial sub-range. The equilibrium
value is given by

σeq =

b
∆
∆

!−4/3

p
−4/3
=
24 c f il

A-POSTERIORI ANALYSISThe open-source code PARIS1 (Aniszewski et al., 2021)
has been employed to solve the unsteady incompressible
Navier-Stokes equations. It uses a second-order Runge-Kutta
technique for time integration and spatial discretization is realized by the finite-volume approach on a regular, cubic staggered grid with second-order centered difference schemes. In
the framework of the projection method, the pressure fieldis calculated by a multi-grid Poisson solver provided by the
HYPRE library.
The Taylor-Green vortex (Brachet et al., 1983) is a challenging test case for laminar-turbulent transition and this configuration consists of a cube with side length 2π and periodic
boundaries in all directions. The velocity field is initialized as(8)u(x, y, z) = cos(x) sin(y) sin(z)(9)

v(x, y, z) = − sin(x) cos(y) sin(z)(10)

w(x, y, z) = 0(11)Referring to the initial state, the Reynolds number is 1600 andthe density is assumed to be constant. Two different resolutions have been investigated, i.e. the benchmark DNS and amuch coarser 323 grid for all a-posteriori LES runs. Hence,the computational costs are considerably different. The nondimensional simulation time ranges from t = 0 to t = 25 andthe corresponding development of the flow in the DNS isshown in Fig. 4. Through vortex breakdown, the flow evolves
from a quasi-laminar initial condition to fairly homogeneous
fully developed turbulence at the final stage considered.
Figure 5 shows the temporal evolution of the mean kinetic energy and its dissipation rate for the no-model LES, theSmagorinsky model, the Bardina/Liu model and the sensorbased mixed model. Both the filtered and unfiltered DNS results are included as references. As expected, the EV model isoverly dissipative during the quasi-laminar initial stage. The
SS model underpredicts the peak dissipation and shows firstb has been
where the ratio of explicit-to-implicit filter width ∆/∆further reduced to the filter coefficients following Lund (1997).
For the explicit test filter applied here, c f il = 1/12 yields σeq =2−2/3 ≈ 0.63.

A-PRIORI ANALYSIS
To justify the above mixed modeling idea, the SGS energy transfer behavior of the separate functional and structuralbase models is analyzed by means of an a-priori analysis. The
comparably well-defined state of homogeneous isotropic turbulence (closely resembled by the final stage of the TaylorGreen vortex, cf. t = 25 in Fig. 4) is chosen for this purpose.
The Direct Numerical Simulation (DNS) database, uniformly
discretized by 5123 grid points, has been explicitly filteredusing a Gaussian filter kernel for varying normalized filter
width ∆/∆DNS . It is particularly instructive to analyze the performance conditional on the coherent structure function FCS
due to the intricate dissipation-enstrophy interplay in turbulent

1 http://www.ida.upmc.fr/∼zaleski/paris/2

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

Figure 1. Homogeneous isotropic turbulence: Means of the sub-grid activity sensor Θ conditional on the coherent structure function
FCS , where −1 and +1 correspond to pure strain and pure rotation, respectively. Results are shown for filter width ∆/∆DNS = 3 (blue),
5 (red), 9 (yellow), 13 (purple), 17 (green). Original formulation based on enstrophy on the left, present formulation based on enstrophyand dissipation on the right.

Figure 2. Homogeneous isotropic turbulence: Probability density function (PDF) of the coherent structure function FCS , where −1
and +1 correspond to pure strain and pure rotation, respectively. Results are shown for filter width ∆/∆DNS = 3 (blue), 5 (red), 9(yellow), 13 (purple), 17 (green).

signs of instability, i.e. oscillations of the dissipation curve,
during the high-dissipation phase around 5 < t < 15. In contrast, the sensor-based mixed model correctly reduces to thefiltered DNS in the initial stage, improves the peak dissipation
and also appears to be more robust during the high-dissipation
phase. A clear improvement compared to the no-model LES
can be observed as well.
The difference between the original sensor by Chapelier
et al. (2018) based on enstrophy (Ωi j Ωi j ) and the present modified sensor based on the sum of enstrophy and dissipation(E = (Ωi j Ωi j + Si j Si j )/2) can be discerned from Fig. 6. Although the overall results are quite similar, the mixed modelusing the newly proposed sensor is somewhat less oscillatory,
which can be seen especially during the high-dissipation phasearound 5 < t < 15. Only the modified sensor variant based on
E is thus used in the following.
A viable alternative for the structural part of the mixed
model is the gradient model by Clark et al. (1979):
τiSS
j =

∆2 ∂ ui ∂ u j
12 ∂ xk ∂ xk

7. As expected, the curves start to diverge considerably when
the level of under-resolution becomes significant at t ≈ 4. In
this case, the improvement through the sensor-based regularization is even more obvious than with the Bardina/Liu model
combination (Fig. 5).

CONCLUSIONS AND OUTLOOKA strategy for mixed modeling in LES is proposed bymeans of a sub-grid activity sensor based blending between thefunctional and structural base models. This mixed model outperforms the separate base models (here: Bardina/Liu or Clarkand Smagorinsky) for the Taylor-Green vortex test case and isoscillation-free without additional regularization like averaging (in homogeneous direction), relaxation in time or clipping(of backscatter). It is worth noting that the overall model isalso parameter-free, apart from the choice of the test filter.
Since not discussed here, future work will focus on thewall treatment for the proposed mixed model. This can beachieved either by replacing the current base models withmodels that already incorporate the correct wall scaling, e.g.
(Nicoud et al., 2011), or by using a wall-scaling sensor, e.g.
(Hasslberger et al., 2021), in addition to the non-wall scaling
sensor Θ for blending. Note that also the structural models by(12)

Again, the sensor-based mixed model leads to a clear improvement over the separate structural model as can be seen in Fig.

3

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

Figure 3. Homogeneous isotropic turbulence: Means of the “true” SGS energy transfer τiDNS
j Si j (EpsTau, continuous lines) and modelenergy transfer τimodelS(EpsMod,markersymbols)
conditionalonthecoherent
structurefunction
FCS . Eps- and Eps+ indicate forwardijjscatter ε − and backward scatter ε + . Results are shown for the Bardina/Liu model (top row) and Smagorinsky model (bottom row) for
moderate, i.e. ∆/∆DNS = 9 (left), and large filter width, i.e. ∆/∆DNS = 33 (right), respectively. Note the different scales of the ordinate
axes.

Figure 4. Stages of the Taylor-Green vortex as seen in the reference DNS: Instantaneous views of FCS = 0 iso-contours coloured by
the velocity magnitude at simulation times t = 2.5, 5, 7.5, 10, 15, 25.

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

0.1

0.01

0.08
0.06

0.008
0.006

0.04

0.004

0.02

0.002

0

0

5

10

15

20

DNS
filtered DNS
no-modelSmagorinsky
Bardina/Liu
mixed,sensor

0.012

dissipation of KE

0.12

kinetic energy KE

0.014

DNS
filtered DNS
no-modelSmagorinsky
Bardina/Liu
mixed,sensor

0

25

0

5

10

t

15

20

25

t

Figure 5. Taylor-Green vortex: Volume-averaged kinetic energy (left) and its dissipation rate (right) versus non-dimensional time forthe reference DNS, the no-model LES, the Smagorinsky model, the Bardina/Liu model and the sensor-based mixed model.
0.14
0.12

0.01

0.08
0.06

0.008
0.006

0.04

0.004

0.02

0.002

0

0

5

10

15

20

DNS
filtered DNS
no-modelmixed,sensor,original
mixed,sensor,present

0.012

dissipation of KE

0.1
kinetic energy KE

0.014

DNS
filtered DNS
no-modelmixed,sensor,original
mixed,sensor,present

0

25

0

5

10

t

15

20

25

t

Figure 6. Taylor-Green vortex: Volume-averaged kinetic energy (left) and its dissipation rate (right) versus non-dimensional time forthe reference DNS, the no-model LES, the mixed model using the original sensor based on enstrophy only and the mixed model using
the newly proposed sensor based on the sum of enstrophy and dissipation.
0.14
0.12

0.01

0.08
0.06

0.008
0.006

0.04

0.004

0.02

0.002

0

0

5

10

15

20

DNS
filtered DNS
no-modelClark
mixed,sensor

0.012

dissipation of KE

0.1
kinetic energy KE

0.014

DNS
filtered DNS
no-modelClark
mixed,sensor

0

25

t

0

5

10

15

20

25

t

Figure 7. Taylor-Green vortex: Volume-averaged kinetic energy (left) and its dissipation rate (right) versus non-dimensional time forthe reference DNS, the no-model LES, the separate model by Clark et al. and the sensor-based mixed model using the structural modelby Clark et al. instead of Bardina/Liu et al.

Bardina/Liu et al. as well as Clark et al. suffer from incorrectnear-wall scaling (Silvis et al., 2017).

R et al. 2021 Parallel, robust, interface simulator (PARIS).
Computer Physics Communications 263, 107849.
Bardina, J 1983 Improved turbulence models based on large
eddy simulation of homogeneous, incompressible, turbulent
flows. PhD thesis, Stanford University.
Brachet, ME, Meiron, DI, Orszag, SA, Nickel, BG, Morf, RH& Frisch, U 1983 Small-scale structure of the Taylor–Greenvortex. Journal of Fluid Mechanics 130, 411–452.
Chapelier, J-B, Wasistho, B & Scalo, C 2018 A coherent vorticity preserving eddy-viscosity correction for large-eddy
simulation. Journal of Computational Physics 359, 164–182.

REFERENCES
Anderson, BW & Domaradzki, JA 2012 A subgrid-scalemodel for large-eddy simulation based on the physics ofinterscale energy transfer in turbulence. Physics of Fluids
24 (6), 065104.
Aniszewski, W, Arrufat, T, Crialesi-Esposito, M, Dabiri, S,
Fuster, D, Ling, Y, Lu, J, Malan, L, Pal, S, Scardovelli,512th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19–22, 2022
larity subgrid-scale models as deduced from measurements
in a turbulent jet. Journal of Fluid Mechanics 275, 83–119.
Lund, TS 1997 On the use of discrete filters for large eddy
simulation. Annual Research Briefs pp. 83–95.
Nicoud, F, Toda, HB, Cabrit, O, Bose, S & Lee, J 2011 Using
singular values to build a subgrid-scale model for large eddy
simulations. Physics of Fluids 23 (8), 085106.
Silvis, MH, Remmerswaal, RA & Verstappen, R 2017 Physical consistency of subgrid-scale models for large-eddy simulation of incompressible turbulent flows. Physics of Fluids
29 (1), 015105.
Smagorinsky, J 1963 General circulation experiments withthe primitive equations: I. the basic experiment. Monthly
Weather Review 91 (3), 99–164.

Clark, RA, Ferziger, JH & Reynolds, WC 1979 Evaluation ofsubgrid-scale models using an accurately simulated turbulent flow. Journal of Fluid Mechanics 91 (1), 1–16.
Hasslberger, J, Engelmann, L, Kempf, A & Klein, M 2021 Robust dynamic adaptation of the Smagorinsky model basedon a sub-grid activity sensor. Physics of Fluids 33 (1),015117.
Klein, M, Ketterl, S, Engelmann, L, Kempf, A & Kobayashi,
H 2020 Regularized, parameter free scale similarity typemodels for large eddy simulation. International Journal of
Heat and Fluid Flow 81, 108496.
Kobayashi, H 2018 Improvement of the sgs model by usinga scale-similarity model based on the analysis of sgs force
and sgs energy transfer. International Journal of Heat andFluid Flow 72, 329–336.
Liu, S, Meneveau, C & Katz, J 1994 On the properties of simi-

6