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

LIFE CYCLE OF TURBULENT KINETIC ENERGY IN RAYLEIGH-BÉNARD
CONVECTION AT MODERATE RA

Myoungkyu Lee
Department of Aerospace Engineering and Mechanics
The University of Alabama
Tuscaloosa, AL 35487, USA
mklee@ua.edu

ABSTRACT

duction mechanism of vertical velocity fluctuations. Yet, theunderstanding of the life-cycle of large-scale motions is still
incomplete. The objective of the current work is to study thelife-cycle of TKE in Rayleigh-Bénard convective flows.
In this work, we investigate the terms in the TKE budget
equation in terms of length-scale and wall-normal locations
to understand the life-cycle of TKE in Rayleigh-Bénard convection problems at different Ras. Among numerous earlier
works in this context, we highlight a few works. (Deardorff &Willis, 1967) performed TKE budget analysis with the flows
up to Ra = 2.5 × 106 . Kerr (2001) studied the effect of aspectratio on the TKE dissipation. More recently, Krug et al. (2020)
performed the spectral analysis on the production of TKE andrevealed the importance of LSMs in Rayleigh-Bénard convection at high Ra. Separately, Lee & Moser (2019) performed
spectral analysis of the budget equations of velocity varianceswith DNS of channel flow at high Re. The same spectral analysis technique is used for this study and summarized in the next
section.

The life-cycle of turbulent kinetic energy (TKE) in the
canonical Rayleigh-Bénard flows is investigated. The spectral densities of terms in the TKE budget equation are computed with direct numerical simulations (DNSs) at Rayleighnumbers, Ra, up to 4 × 108 . The result exhibits that both production and dissipations are most significant at the center between two planes. Small-scale motions are dominant for dissipation, while the non-negligible amount of TKE is dissipatedby large-scale motions in the near-wall region. On the other
hand, large-scale motions are dominant for the production ofTKE, but small-scale motions also contribute to the production of TKE. It suggests that the traditional energy cascade
model, which describes the energy transfer from large-scale
to small-scale, has a significant role in the life-cycle of TKE.
Furthermore, the cascade process from large-scale motion tosmall-scale motion involves transports by nonlinear mechanisms, such as turbulence and pressure-velocity correlations.
The length-scale of the wall-normal transport mechanism inthe overlap region grows linearly with wall-normal distances,and such linear growth regions increase with Ra.

METHOD
DNSs of the incompressible Rayleigh-Bénard convection
are performed for this study. We solve the normalized governing equations for the incompressible Rayleigh-Bénard convection. Those are

INTRODUCTION
The super-structure or very large-scale motions (VLSMs)of turbulent flows are observed in various types of wallbounded turbulence. Such large-scale motions are stronger
in boundary-driven wall-bounded turbulence than in pressuregradient-driven wall-bounded turbulences (Lee & Moser,
2015; Pandey et al., 2018; Lee & Moser, 2018). RayleighBénard convective flow is another boundary condition driven
flow with VLSMs (Stevens et al., 2018; Krug et al., 2020;
Green et al., 2020; Blass et al., 2021). The large-scale motions(LSMs) of Rayleigh-Bénard convective flows are important invarious engineering and scientific problems, e.g., atmosphericflows, boiling devices, etc., since they govern the macroscopictransport of mass and turbulent kinetic energies (Ahlers et al.,2009).
Recently, some properties of the large-scale motions inRayleigh-Bénard convective flows have been revealed by several DNS studies. The TKE spectral densities at different
Ra have peaks at k = 1 where k is the circular wavenumber(Stevens et al., 2018), and the size of large-scale motions growwith Ra (Krug et al., 2020). Krug et al. also reported thatthe size of large-scale structures in the temperature field is different from the size of large-scale motions in vertical velocity fields, and the difference is strongly connected to the pro-

∂ uk
=0
∂ xk
∂ ui
∂ ui uk
=−
+
∂t
∂ xk

r

Pr ∂ 2 ui∂p−+ θ δi3Ra ∂ xk ∂ xk ∂ xi

(1)∂θ
∂ θ uk
∂ 2θ
1
=−
+√
∂t
∂ xk
PrRa ∂ xk ∂ xk
where Pr and Ra are Prandtl number and Rayleigh number, respectively. Also, ui , p, θ and δi j are the velocity components,pressure, normalized temperature and Kronecker’s delta, respectively. x1 and x2 are homogeneous directions, and x3 isthe wall-normal direction. No-slip and no-penetration conditions for velocity fields are applied at top and bottom walls.
Boundary conditions for normalized temperatures are1θx3 =0 = ,21

θx3 =h = −

1
2

(2)

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

Ra

L/H

N∥

N⊥

Nu

Re∥

Re⊥

ReT

δθ

1E5

1 × 105

10π

512

96

4.34

55.6

40.1

68.6

1.152 × 10−1

6E6

6 × 106

5π

768

256

13.62

425.3

299.4

520.1

3.672 × 10−2

4E8

4 × 108

5π

2048

512

46.37

3141.7

1962.43704.3

1.078 × 10−2

Table 1. Simulation parameters. L: simulation domain size in the homogeneous directions, H: distance between walls, N∥ : numberof grid points in spectral domain, N⊥ : number of collocation points in wall-normal direction, Re∥ : horizontal Reynolds number, Re⊥ :wall-normal Reynolds number, ReT : total Reynolds number, δθ : thermal boundary layer thickness∂ ⟨ũ′i u′i u′1 ⟩ ∂ ⟨u′i ũ′i ũ′1 ⟩ ∂ ⟨ũ′i u′i u′2 ⟩ ∂ ⟨u′i ũ′i ũ′2 ⟩
−
+
−
∂ r1
∂ r1
∂ r2
∂ r2


1 ∂ ⟨ũ′i u′i u′3 ⟩ ∂ ⟨u′i ũ′i ũ′3 ⟩
−
+
2
∂ x3
∂ x3


′ 
∂
u
∂ ũ′
+ ũ′i ũ′3 i + u′i u′3 i
∂ x3
∂ x3
r
!
 ′
∂ ui ∂ ũ′iPr ∂ ⟨u′i ũ′i ⟩ ∂ ⟨u′i ũ′i ⟩
2Rε =
+
+
Ra
∂ x3 ∂ x3
∂ r12
∂ r22
∥

where h is the distance between two walls.
We usethe velocity-vorticity formulation from Kim et al. (1987).
Spectral-Galerkin method and seventh-order basis splinecollocation method is used for derivatives in homogeneous andthe wall-normal directions, respectively. Orszag’s 3/2 rules areapplied with zero padding to remove the aliasing for computing quadratic products (Orszag, 1971). The low-storage thirdorder Runge-Kutta method is used for the time integration.
Modified in-house simulation code, which has been used for
Lee & Moser (2015, 2018), is used for simulations. Consult
Lee et al. (2013, 2014) for details of the simulation code implementations.

2RT =The Fourier transforms of terms in (4) with respect to theseparation distances provide spectral densities of terms in the
TKE budget equation.

We simulate three cases at different Ras with fixed Pr,
Pr = 1. The Nusselt number, Nu and thermal boundary layer
thickness, δθ , of three cases are in table 1. The values in table 1 agrees values from Krug et al. (2020) within statisticaluncertainties.

∂ EK∥= EP + ET + ET⊥ + EΠ + ED − Eε∂t

The two-point correlation function corresponding to theturbulent kinetic energy (TKE, K = ⟨ui ui ⟩/2) is defined withquantities separated in homogeneous directions, r, e.g. ũi =ui (x + r) where x = (x1 , x2 , x3 ) and r = (r1 , r2 , 0).

2RK (r1 , r2 , x3 ) = ⟨u′i ũ′i ⟩The integration of terms in (5) over k1 and k2 are the terms in
the TKE budget equations (Deardorff & Willis, 1967).
∂K= P+T +Π+D−ε∂t(6)where P, T , Π, D and ε are production, wall-normal transport
by turbulence, wall-normal transport by pressure-velocity interaction, wall-normal transport by viscosity and dissipation,∥
respectively. Note that the integration of ET over k1 and k2 is
zero.

(3)

where ⟨·⟩ denotes the averaging in time and homogeneous directions. With some mathematical manipulations, the budget
equation of RK is obtained as follows.

∥

lim RT =

r1 →0
r2 →0

∂ RK
∥
= RP + RT + RT⊥ + RΠ + RD − Rε∂t(5)

Z ∞ Z ∞−∞ −∞∥ET dk1 dk2 = 0,

∀x3

(7)

∥

(4)

Therefore, we interpret ET as the term for inter-scale transfer
at any fixed wall-normal distance, and this term is not shown in
the TKE budget equation, eq (6). Also, both T and Π are wallnormal transport terms with nonlinear interaction of velocity
components, while D is the linear wall-normal transport terms.
Note that the pressure field in incompressible flows is the resultof nonlinear interactions of velocity components to satisfy thedivergence-free condition. Therefore, we combine T and Π asfollows.

where

2RP =⟨u′3 θ̃ ′ ⟩ + ⟨ũ′3 θ ′ ⟩1 ∂ ⟨ũ′i u′i u′3 ⟩ ∂ ⟨u′i ũ′i ũ′3 ⟩2RT⊥ = −
+
2
∂ x3
∂ x3
∂ ⟨u′3 p̃′ ⟩ ∂ ⟨ũ′3 p′ ⟩
2RΠ =
+
∂ x3
∂ x3
r
2
′
′
Pr ∂ ⟨ui ũi ⟩2RD =Ra ∂ x32

N = T +Π(8)

EN = ET⊥ + EΠ

(9)

and

2

0.030
0.025
0.020
0.015
0.010
0.005
0.000 1
10

100
x3/

K

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

Ra = 1 × 105
Ra = 6 × 106
Ra = 4 × 108
100

x3/

101

10 1

102

0

Z h/2
0

ED dx3 = 0,∀(k1 , k2 )x3/

Note that integration of EN or ED over x3 at any wavenumber
is zero.

EN dx3 =

101
/

100

0.01080.00960.00840.00720.00600.00480.00360.00240.00120.0000

(10)

100

10 1

(b) Ra = 6 × 106
101

100

Finally,

102

/

101

(11)

x3/

∂ EK
∥
= EP + ET + EN + ED − Eε
∂t

102

101

Figure 1. One-dimensional profile of turbulent kinetic energy (TKE)Z h/2(a) Ra = 1 × 105

RESULT
In this section, the profile and spectral density of TKE are
discussed. Then, its budget terms are discussed in the order of
production, dissipation, linear wall-normal transport, nonlinear wall-normal transport, and inter-scale transfer.
The one-dimensional profile of TKE is shown in the figure 1. Locations of TKE peaks increase with Ra, while the values of TKE peaks decrease with Ra. Also, the increasing rateof TKE at the near-wall region exhibit some Ra dependencies.
Corresponding one-dimensional spectral densities of TKE are
shown in the figure 2. Note that figures with spectral densities
are represented as functions of k and x3 /δθ . k is premultiplied
to all spectral densities to compensate for the log scales in k
domains. The premultiplied TKE spectra shows strong peakat x3 /δθ = 1 and λ /δθ = 40 (λ /h ≈ 4.5) at Ra = 1 × 105 .
Interestingly, LSMs at Ra = 6 × 106 show bimodal peaks, andthe peaks show scale separations as Ra increases. It is possible that this observation is an artifact of arbitrary bin size toobtain the spectral density in k with ∆k = 0.5. There is only
one more data point between two peaks at Ra = 6 × 106 , butthere are two data points between two peaks at Ra = 4 × 108 .
Nonetheless, the data confirms that the length-scale of superstructure rolls grows with Ra, which is consistent with (Pandey
et al., 2018). Additional simulation results with much bigger
simulation domains at higher Ra are necessary to verify thisobservation, but the costs of such simulations are prohibitively
expensive.
One-dimensional profiles of production and dissipation of
TKE in figure 3 show that most of the produced TKE get dissipated at the same wall-normal locations, and a small portion of
TKE which is not dissipated in the outer region is dissipated at
the near-wall region, x3 /δθ < 1. Here, x3 /δθ is also premultiplied to the budget terms to compensate the log-scale in x3 /δθsimilar to k-premultiplied spectral density of TKE. The imbalance between production and dissipation is more vivid in the

100
10 1

(c) Ra = 4 × 108
101

/

102

0.027
0.024
0.021
0.018
0.015
0.012
0.009
0.006
0.003
0.000

103

0.00810.00720.00630.00540.00450.00360.00270.00180.00090.0000(x3/ )B

Figure 2. k-Premultiplied spectral densities of TKE

0.12
0.10
0.08
0.06
0.04
0.02
0.00 110

Ra = 1 × 105
Ra = 6 × 106
Ra = 4 × 108
PK
K

100

x3/

101

102

Figure 3. One-dimensional profile of TKE production (PK )and dissipation (εK )spectral density of those. The spectral densities of production
have peaks at the middle point between walls, and corresponding length-scales of peaks increase with Ra even with outerscaling, λ /h. Similarly, the spectral densities of dissipation
have peaks at the middle point between walls. However, thelength-scale of dominant dissipation is λ /δθ = 7 − 8 and it is3

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

0.06480.05760.05040.04320.03600.02880.02160.01440.00720.0000(a) Ra = 1 × 105
101
/

102

0.036
0.032
0.028
0.024
0.020
0.016
0.012
0.008
0.004
0.000

x3/

101

100

10 1

(b) Ra = 6 × 106
100

101

102

/

x3/

101
100
10 1

(c) Ra = 4 × 108
101

/

102

103

(a) Ra = 1 × 105
100

101
/

102

0.03240.02880.02520.02160.01800.01440.01080.00720.00360.0000

101

x3/

100

10 1

100

10 1

0.0486
0.04320.03780.03240.02700.02160.01620.01080.00540.0000(b) Ra = 6 × 106
100

101

102

/

101
x3/

10 1

100
x3/

x3/

100

0.020250.01800
0.01575
0.01350
0.01125
0.00900
0.00675
0.00450
0.002250.00000

100
10 1(c) Ra = 4 × 108
101

/

102

103

0.04050.03600.03150.02700.02250.01800.01350.00900.00450.0000Figure 5. k- and x3 -Premultiplied spectral density of εK

Figure 4. k- and x3 -Premultiplied spectral density of PK

almost invariant to Ra. Interestingly, a small amount is dissipated at the near-wall region, and their length-scale is consistent with the large-scale structure of production. The different spectral behaviors of production and dissipation imply thatthe roles of wall-normal transport mechanisms and inter-scaletransfer are important in the life-cycle of TKE.
Figure 6 shows that the non-linear wall-normal mechanism transports TKE from the outer region to the near-wall region. Then, the linear wall-normal mechanism transports TKE
further to the wall. However, the spectral behaviors of wallnormal transport mechanisms are more complex. The observation with one-dimensional profiles still holds with large-scale
motions. (Figure 7) When the length scale gets smaller, thenonlinear wall-normal mechanism also transports some TKE
away from the wall. Also, the length scale of the donor region in the nonlinear mechanism linearly increases with wallnormal distance, λ ∼ x3 , and this linear growth is vivid when
Ra increases. There is a small region where the linear mechanism also transports TKE away from the wall (Figure 8). Thisbehavior is only shown in the region x3 /δθ > 1, and the magnitude is negligible. Generally, the linear mechanism transports
TKE, which is transported from the outer region by the nonlinear mechanism, to the near-wall region in all length scales.
The spectral behavior of inter-scale TKE transfer is shownin figure 9. At Ra = 1 × 105 , the strong donor region is where
the production spectra show the peak. A small inverse energy
cascade is observed in the near-wall region at Ra = 1 × 105 ,but the magnitude of the inverse TKE cascade is small. Also,
the inverse energy cascade is not observed at higher Ra flows.
The inter-scale transfer mechanism in the outer region, say
x3 /h > 0.1, shows different behavior from the near-wall region. The donor region has only one peak at Ra = 1 × 105 .
However, two peaks are observed in the donor region at outerflow in higher Ra flows, while the recipient region has a single
peak in all flows. Among the peaks in the donor region at highRa flows, the peak with larger-scale is at the same length-scale
where the production mechanism has peaked. However, thepeaks at the intermediate length scale are at different lengthscales where the production spectra have peaks. Instead, theTKE transferred from the intermediate length scale is from theregion closer to the wall. The TKE produced at the large-scale
motion at the middle point between walls transported toward
walls, say x3 /h < 0.1 by the nonlinear wall-normal transportmechanism. Then, the inter-scale mechanism transfers the energy from the large-scale to the intermediate length scale. Finally, the nonlinear wall-normal mechanism transports TKE
toward the middle point between walls at a given length scale.

4

0.03
0.02
0.01
0.00
0.01
0.02
0.03 110

Ra = 1 × 105
Ra = 6 × 106
Ra = 4 × 108

K
DK

100
x3/

(x3/ )B

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

100

101

x3/

10 1

102

x3/
100

101

x3/

101

0.00330.00000.00330.00660.00990.01320.01650.01980.0231100

10 1(b) Ra = 6 × 106
100

101

x3/

101
100
10 1

(c) Ra = 4 × 108
101

/

102

103

0.00150.00100.00050.00000.00050.00100.00150.00200.00250.0030(b) Ra = 6 × 106
100

101

102

/

0.000990.00066
0.00033
0.00000
0.00033
0.00066
0.000990.001320.00165100
10 1(c) Ra = 4 × 108
101

/

102

103

Figure 8. k- and x3 -Premultiplied spectral density of DK

102

/

102

101

102

/

101
/

100

10 1

x3/

10 1

x3/

0.01280.00640.00000.00640.01280.01920.02560.03200.03840.0448(a) Ra = 1 × 105

100

101

Figure 6. One-dimensional profile of nonlinear wall-normaltransport of TKE (NK ), and linear wall-normal transport ofTKE (DK ).

100(a) Ra = 1 × 105

0.004
0.002
0.000
0.002
0.004
0.006
0.008
0.010
0.012
0.014

cess of TKE from the large-scale to the small-scale, some TKEgets transported toward walls and then transported away fromthe wall.
In the early part of this work, we speculate the existenceof bi-modal spectral peaks of TKE at high Ra flows, but it isunable to verify it because it could be an artifact of a small simulation domain. Even though it is at the different wall-normal
locations, the inter-scale transfer also shows the bimodal peaksin the donor region at high Ra flows. We speculate that these
two bimodal behavior at high Ra flows is connected, but, again,
we may need to have simulation results with a larger simulation domain and perhaps higher Ra flows.
There are several future works. First, we need a better
scaling parameter than δθ and h. δθ successfully describes thescaling of dissipation at the middle point between walls. Other
than that, we are unable to find any invariant scaling with δθ
in flows at different Ras even in the near-wall region. Second,
this work is limited to the TKE, but each component of TKE
should be studied. Especially, the wall-normal velocity component is the only one with a non-zero production mechanism.
Therefore, inter-component energy transfer should have an important role in the life-cycle of TKE. Third, we can apply thesame analysis techniques to the temperature fluctuations andvelocity-temperature correlations. The velocity-temperaturecorrelation is particularly important because it is a production

0.01100.00550.00000.00550.01100.01650.02200.02750.0330Figure 7. k- and x3 -Premultiplied spectral density of NK

CONCLUSION
In this work, we apply spectral analysis to the terms in the
TKE transport equation. The production and dissipation aremostly concentrated at the middle point between walls. However, their dominant length scales are different. The inter-scaletransfer mechanism is mostly responsible for such scale separation between production and dissipation. In the cascade pro-5

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

0.002250.00000
0.002250.00450
0.00675
0.00900
0.01125
0.01350
0.01575
0.01800

x3/

100

10 1

(a) Ra = 1 × 105
100

101
/

102

0.01800.01440.01080.00720.00360.00000.00360.00720.01080.0144

x3/

101

100

10 1

(b) Ra = 6 × 106
100

101

102

/

x3/

101
100
10 1

(c) Ra = 4 × 108
101

/

102

Rayleigh-Bénard convection. Reviews of Modern Physics
81 (2), 503–537.
Blass, Alexander, Verzicco, Roberto, Lohse, Detlef, Stevens,
Richard J.A.M. & Krug, Dominik 2021 Flow organizationin laterally unconfined Rayleigh-Bénard turbulence. Journal of Fluid Mechanics 906, A26.
Deardorff, J. W. & Willis, G. E. 1967 Investigation of turbulent thermal convection between horizontal plates. Journalof Fluid Mechanics 28 (4), 675–704.
Green, Gerrit, Vlaykov, Dimitar G, Mellado, Juan Pedro &Wilczek, Michael 2020 Resolved energy budget of superstructures in Rayleigh–Bénard convection. Journal of Fluid
Mechanics 887, A21.
Kerr, R. M. 2001 Energy Budget in Rayleigh-Bénard Convection. Physical Review Letters 87 (24), 244502.
Kim, John, Moin, Parviz & Moser, Robert 1987 Turbulence
statistics in fully developed channel flow at low reynolds
number. Journal of Fluid Mechanics 177, 133–166.
Krug, Dominik, Lohse, Detlef & Stevens, Richard J. A. M.
2020 Coherence of temperature and velocity superstructures in turbulent rayleigh-bénard flow. Journal of Fluid
Mechanics 887, A2.
Lee, Myoungkyu, Malaya, Nicholas & Moser, Robert D. 2013
Petascale direct numerical simulation of turbulent channel
flow on up to 786k cores. In SC13, the International Conference for High Performance Computing, Networking, Storage and Analysis, pp. 1–11.
Lee, Myoungkyu & Moser, Robert D. 2015 Direct numerical simulation of turbulent channel flow up to Reτ ≈ 5200.
Journal of Fluid Mechanics 774, 395–415.
Lee, Myoungkyu & Moser, Robert D. 2018 Extreme-scale motions in turbulent plane couette flows. Journal of Fluid Mechanics 842, 128–145.
Lee, Myoungkyu & Moser, Robert D. 2019 Spectral analysisof the budget equation in turbulent channel flows at highreynolds number. Journal of Fluid Mechanics 860, 886–
938.
Lee, Myoungkyu, Ulerich, Rhys, Malaya, Nicholas & Moser,
Robert D. 2014 Experiences from leadership computing insimulations of turbulent fluid flows. Computing in Science
Engineering 16 (5), 24–31.
Orszag, Steven A. 1971 On the elimination of aliasing in finitedifference schemes by filtering high-wavenumber components. Journal of the Atmospheric Sciences 28, 1074.
Pandey, Ambrish, Scheel, Janet D. & Schumacher, Jörg 2018
Turbulent superstructures in rayleigh-bénard convection.
Nature Communications 9 (1), 1–11.
Stevens, Richard J A M, Blass, Alexander, Zhu, Xiaojue, Verzicco, Roberto & Lohse, Detlef 2018 Turbulent thermal superstructures in Rayleigh-Bénard convection. Physical Review Fluids 3 (4), 041501(R).

103

0.03600.02880.02160.01440.00720.00000.00720.01440.02160.0288∥Figure 9. k- and x3 -Premultiplied spectral density of TK

mechanism of TKE. Yet, the analysis of velocity-temperaturecorrelation is challenging because it requires studies not only
on the magnitude but also on the phase in the spectral densitiesof budget terms. Finally, we limited our simulation cases withunity Pr in this work. The effect of Pr must be investigated.

ACKNOWLEDGEMENT
This research used resources of the Argonne Leadership
Computing Facility, which is a DOE Office of Science User
Facility supported under Contract DE-AC02-06CH11357.

REFERENCES
Ahlers, Guenter, Grossmann, Siegfried & Lohse, Detlef
2009 Heat transfer and large scale dynamics in turbulent

6