New J. Phys. 19 (2017)065006 https://doi.org/10.1088/1367-2630/aa6fe8
PAPER
Extreme-value statistics from Lagrangian convex hull analysis for
homogeneous turbulent Boussinesq convection and MHD
convection
J Pratt
1
, A Busse
2
, W-C Müller
3
, N W Watkins
4,5,6,7
and S C Chapman
4,8,9
1
Astrophysics Group, University of Exeter, Exeter EX4 4QL, United Kingdom
2
School of Engineering, University of Glasgow, Glasgow G12 8QQ, United Kingdom
3
Center for Astronomy and Astrophysics, ER 3-2, TU Berlin, Hardenbergstr. 36, D-10623 Berlin, Germany
4
Centre for Fusion, Space and Astrophysics, Physics Department, University of Warwick, Coventry CV4 7AL, United Kingdom
5
Institut für Physik und Astronomie, Universität Potsdam, Campus Golm, Haus 28, Karl-Liebknecht-Strasse 24/25,
D-14476 Potsdam-Golm, Germany
6
Centre for the Analysis of Time Series, London School of Economics and Political Science, London, United Kingdom
7
Faculty of Science, Technology, Engineering and Mathematics, Open University, Milton Keynes, United Kingdom
8
Max-Planck-Institut für Physik komplexer Systeme, D-01187 Dresden, Germany
9
Department of Mathematics and Statistics, University of Tromso, Norway
E-mail: [email protected]
Keywords: turbulence, magnetohydrodynamics, Lagrangian statistics, magnetoconvection, turbulent transport
Abstract
We investigate the utility of the convex hull of many Lagrangian tracers to analyze transport properties
of turbulent flows withdifferent anisotropy. In directnumerical simulations ofstatistically
homogeneous and stationary Navier–Stokes turbulence, neutral fluid Boussinesq convection, and
MHD Boussinesq convection a comparison with Lagrangian pair dispersion shows that convex hull
statistics capture the asymptotic dispersive behavior of a large group of passive tracer particles.
Moreover, convex hull analysis provides additional information on the sub-ensemble of tracers that
on average disperse most efficiently in the form of extreme value statistics and flow anisotropy via the
geometric properties of the convex hulls. We use the convex hull surface geometry to examine the
anisotropy that occurs in turbulent convection. Applying extreme value theory, we show that the
maximal square extensions of convex hull vertices are well described by a classic extreme value
distribution, the Gumbel distribution. During turbulent convection, intermittent convective plumes
grow and accelerate the dispersion of Lagrangian tracers. Convex hullanalysis yields information that
supplements standard Lagrangian analysis of coherent turbulent structures and their influence on the
global statistics of the flow.
1. Introduction
Turbulent transport governs the spreading of contaminants in the environment, mixing of chemical
constituents in combustion engines or in stellar interiors, accretion in proto-stellar molecular clouds,
acceleration of cosmic rays, and escape of hot particles from fusion machines. Because of its wide relevance, a
fundamental characterization of the dispersive properties of turbulent flows is of practical interest to physicists
and engineers. Here we examine the broadly relevant case of dispersion of Lagrangian tracer particles in
statistically homogeneous but not necessarily isotropic turbulence.
The Lagrangian viewpoint is particularly suited to the investigation of transport in turbulent fluids. A
Lagrangian description of turbulence is based on following the paths of passive tracer particles in a turbulent
flow. Single-particle diffusion, as originally addressed by Taylor [1], provides a basic characterization of a flow’s
transport properties [2]. A more complete characterization of the turbulent transport has conventionally been
formed from the relative dispersion of two, three, or four particles [3–11]. However, in astrophysical
OPEN ACCESS
RECEIVED
30 September 2016
REVISED
28 March 2017
ACCEPTED FOR PUBLICATION
27 April 2017
PUBLISHED
20 June 2017
Original content from this
work may be used under
the terms of the Creative
Commons Attribution 3.0
licence.
Any further distribution of
this work must maintain
attribution to the
author(s)and the title of
the work, journal citation
and DOI.
© 2017 IOP Publishing Ltd and Deutsche Physikalische Gesellschaft
environments where the effects of magnetic fields, rotation, or gravity are often significant, the more complex
nature of statistically anisotropic or even inhomogenous nonlinear dynamics warrants additional examination.
Dispersion in dynamically anisotropic systems such as vigorously convecting flows [12–16]where preferred
directions exist and spatially coherent, persistent structures like convective plumes can form, motivates the
present consideration of a complementary diagnostic based on a different Lagrangian concept: the convex hull
[17]of a n-particle group (
n
4).
The convex hull is the smallest convex polygon that encloses a group of particles; two dimensional convex
hulls are pictured in figure 1. Convex hull analysis of turbulent dispersion is similar in spirit to following a drop
of dye as it spreads in a fluid, or following a puff of smoke as it spreads in the air, both classical fluid dynamics
problems [18–20]. A large group of tracer particles can be marked, similarly to adding a drop of dye to a fluid
flow, so that the same particles can be identified at all later times. Using the convex hull, a size for the group of
tracer particles that were marked can be calculated at each time.
The convex hull yields statistical information about a class of Lagrangian particles that is not equivalent to
pre-selected tracer particle groups like particle pairs or tetrads. These standard Lagrangian multi-particle
statistics represent a fixed and unique structural relationship between specific tracer particles. The evolution of
particle-pair structures, expressed e.g. as separation and orientation, is analyzed as the pair of particles is
advected by the fluid. In contrast, the convex hull does not establish a unique link between the tracers that
generate it, but continuously selects from a predefined group, based on which tracer particles have ventured
furthest from the geometrical center of the ensemble. Unlike particle pairs or tetrads, the particles that constitute
the convex hull are dynamically changing. The definition of the convex hull thus corresponds to a filtering based
on the entire dynamical past of each particle in the group. The convex hull captures the extremes of the
excursions of a group of particles, information relevant to the non-Gaussian aspects of the dynamics. The
behavior of particles that do not exhibit the fastest dispersion is filtered by the convex hull, allowing a
classification of particle dynamics with regard to their dispersion efficiency. In this work we begin to explore this
link to extreme value theory, which has the potential to provide new physical insight for turbulent diffusion. The
dynamical relation between the Lagrangian particle population forming the convex hull and the bulk ensemble
of tracer particles enclosed by it represents another aspect of this diagnostic that could be exploited in
investigations of turbulent structure formation.
Figure 1. An illustration of a two dimensional convex hull (solid line)surrounding a group of particles (solid points)as they disperse in
time. The time progression is indicated by arrows, and the particles in each of the three convex hulls shown are the same.
2
New J. Phys. 19 (2017)065006 J Pratt et al
In recent years, convex hull calculations have beenused to study diverse topics suchas the size of spreading
GPS-enabled drifters moving on the surfaces of lakes and rivers [21,22], star-forming clusters [23], forest fires
[24], proteins [25,26], or clusters of contaminant particles [27]. Studies of the relationships between random
walks, anomalous diffusion, extreme statistics and convex hulls have been motivated by animal home ranges
[28–34]. Convex hulls have also been used to study analytical statistics of Burgers turbulence by analogy with
Brownian motion [35–37].
MHD turbulence [38,39]and turbulence during hydrodynamic convection [40–42]are areas where
statistical analysis of Lagrangian particles has begun to be applied only recently. This work presents new
Lagrangian results from three-dimensional direct numerical simulations of turbulent MHD Boussinesq
convection, and compares them with turbulent hydrodynamic Boussinesq convection and homogeneous
isotropic turbulence. It is structured as follows. In section 2we describe the fluid simulations. In sections 3we
present standard Lagrangian pair dispersion and discuss the results of these widely-used statistical tools for
convective flows. In section 4we describe the convex hull analysis that we perform on groups of many
Lagrangian tracer particles. We perform several basic checks on our convex hull calculations. We then compare
the dispersion curves obtained from convex hulls of large groups of particles with the expected scalings for
particle-pair dispersion. In section 5we demonstrate how the convex hull can be used to examine anisotropy. In
section 6we apply extreme value theory, and show that the maximal square extensions of convex hull vertices are
well described by a classic extreme value distribution, the Gumbel distribution. In section 7we summarize the
results of this validation study, and our extreme value statistics. We discuss the potential uses and benefits of
convex hull analysis.
2. Simulations
We investigate three different types of turbulent systems: forced homogeneous isotropic Navier–Stokes
turbulence (simulation NST)[43], Boussinesq convection in a neutral fluid (simulation HC), and Boussinesq
convection in an electrically conducting fluid (simulation MC)[12,44]. These simulations are not designed for
close comparison, but produced for a broad exploration of the convex hull analysis. In each of these direct
numerical simulations, the equations are solved using a pseudospectral method in a cubic simulation volume
with a side of length 2p. The non-dimensional Boussinesq equations for MHD convection in Alfvénic units are
vjB g
t,1
20
wwwnq
¶
¶-´ ´ +´ = -´()ˆ ()
BvB B
t,2
2
h
¶
¶-´ ´ = ()ˆ ()
vv
tT,3
20
qqkq
¶
¶+=-(· ) ˆ(· ) ()
vB0. 4= =·· ()
These equations include the solenoidal velocity field
v
, vorticity
vw=´
, magnetic field
B
, and current
j
B=´ . The quantity θdenotes the temperature fluctuation about a linear mean temperature profile
Tz
0
()
where zis the direction of gravity. In equation (3)this mean temperature gradient provides the convective drive
of the system. In equation (1), the term including the temperature fluctuation θis the buoyancy force. The vector
g
0is a unit vector in the direction of gravity. Three dimensionless parameters appear in the equations:
n
ˆ,
h
ˆ
, and
k
ˆ
. They derive from the kinematic viscosity ν, the magnetic diffusivity η, and thermal diffusivity κ.
For simulation HC, the magnetic field
B
is set to zero. For simulation NST, both magnetic field terms and
temperature terms are zero. A fixed time step and a trapezoidal leapfrog method [45]are used for the time-
integration for simulation NST. The Boussinesq convection simulations HC and MC are integrated in time
using a low-storage 3rd-order Runge Kutta scheme [46]and an adaptive time step, which allows for better time
resolution of large fluctuations that occur during convection.
In this work we discuss turbulent dispersion inan incompressible fluid, where conservation of volume is a
primitive concept. A volume of fluid that is convex at an initial time will occupy the same volume after a period
of dynamic development but will generally change its shape and lose its convexity. Lagrangian tracer particles
that are contained in the initial volume are marked so that they can be followed for the entire time of the
simulation. At any later time, the volume of the convex hull of that group of marked particles is generally not
conserved. This is illustrated in figure 1for a group of particles, and for snapshots taken at three times. The
growth of surface area and volume are natural concepts for convex hulls.
A summary of the fundamental parameters that describe each simulation is given in table 1. In this table, we
define the Reynolds number to be
Re E L
v
12
n
=á ñ ˆ
, where
v2
2
E
v
=
is the kinetic energy, and the brackets
indicate a time-average. We define the characteristic length scale
L
based on the largest-scale motions of the
3
New J. Phys. 19 (2017)065006 J Pratt et al
Table 1. Simulation parameters: grid size N
3
, total number of particles in the simulation
n
p
(10
6
), Reynolds number
R
e
, magnetic Reynolds number
R
em, Prandtl number
P
r
, magnetic Prandtl number
P
rm, Rayleigh number
R
a
,
Kolmogorov microscale kol
h
, Kolmogorov time-scale
t
h
, Lagrangian crossing time
L
CT
, average kinetic energy dissipation rate
v
, Alfvén ratio rA, average Bolgiano–Obukhov length divided by the height of the simulation volume
bo
ℓ
¯,
average number of particles per convex hull
n
pch
, number of convex hulls
N
hulls, initial length scale of convex hull
hull
ℓ
.
N
3
n
p
(10
6
)
R
e
R
em
P
r
P
rm
R
a
(10
5
)kol
h
(10
−3
)
t
h
(10
−2
)
L
CT
(
t
h
)
v
rA
bo
ℓ
¯
n
pch
N
hulls
hull kol
hℓ()
NST 1024
3
3.2 2900 ——— — 4.58 5.25 276 0.15 —— 24 5000 27
HC 512
3
1.0 1500 —2.0 —5 12.6 3.97 340 2.54 —0.28 214 2500 22
MC 512
3
1.0 5100 7650 2.0 1.5 2.22 8.9 2.60 530 4.40 1.78 0.12 48 2500 30
4
New J. Phys. 19 (2017)065006 J Pratt et al
system in question. For statistically homogeneous turbulent convection the characteristic length scale is the
instantaneous temperature gradient length scale TT
0
L
*
=
where T
*
is the root-mean-square (rms)of
temperature fluctuations and ∇T
0
is the constant vertical mean temperature gradient [47]. For non-convective
statistically homogeneous turbulent flows, the characteristic length scale is a dimensional estimate of the size of
the largest eddies,
32
L
Ev
v
=
, where
vk
k22
v
n=á
åñˆ
is the time-averaged rate of kinetic energy dissipation.
The magnetic Reynolds number is defined from the Reynolds number and the magnetic Prandtl number, i.e.
Re Pr Re
mm
=. We measure length in units of the Kolmogorov microscale
kol 3v14
h
n=(ˆ )
and also make use
of the Kolmogorov time-scale
v12
tn=
h
(ˆ )
. The Kolmogorov microscale multiplied by
kmax
, the highest
wavenumber in the simulation, is often used to test whether a simulation is adequately resolved on small spatial
scales. In this work all of the simulations fulfill the standard criterion based on the Kolmogorov microscale
(.k1
5
kol
max h>)for adequate spatial resolution [48]. The Reynolds numbers in table 1are on the order of 10
3
;
the Reynolds numbers and Kolmogorov microscales in table 1are in the same range as current studies of
moderately turbulent flows (e.g. [49,50]).
Formulation of boundary conditions for simulations of turbulent flows is delicate because boundaries
strongly influence the structure and dynamics of the flow. For homogeneous isotropic turbulence, it is standard
to employ boundary conditions that are periodic in x,y, and z. These fully periodic boundary conditions are used
for simulation NST. For convection simulations the choice of fully periodic boundary conditions (also called
homogeneous Rayleigh–Bénard boundary conditions)allows macroscopic elevator instabilities to form [51].
These instabilities destroy the natural pattern of the original turbulent flow field. The convection simulations
discussed in this work use quasi-periodic rather than fully periodic boundary conditions. In quasi-periodic
boundary conditions the only additional constraint is the explicit suppression of mean flows parallel to gravity,
which are removed at each time step. Because our simulations are pseudospectral, the mean flow is
straightforwardly isolated as the zcomponent of the
k
0, 0, 0=(
)
mode in Fourier space, which corresponds to
the volume-averaged velocity in the z-direction. Quasi-periodic boundary conditions combine the conceptual
simplicity of statistical homogeneity with a physically natural convective driving of the turbulent flow. These
boundary conditions do not enforce a large-scale structuring of the turbulent flow, such as the convection-cell
pattern observed when Rayleigh–Bénard boundary conditions are used. In the quasi-periodic simulations
presented in this work, we find no evidence of the macroscopic elevator instability although we follow the
evolution of the flow for long times. Quasi-periodic boundary conditions allow for direct comparison with
simulations that use fully periodic boundary conditions.
In simulation NST the modes
k2.5 3.5<<
are forced using Ornstein–Uhlenbeck processes with a finite
time-correlation on the order of the autocorrelation time of the velocity field (for further details of this forcing
method, see Eswaran and Pope [52]). The convection simulations HC and MC are Boussinesq systems driven
solely by a constant temperature gradient in the vertical direction. The magnetic field present in simulation MC
is generated self-consistently by the flow from a small random seed field through small-scale dynamo action. The
system is evolved until a statistically stationary state is reached. For Boussinesq convection, a length scale that
characterizes the scale-dependent importance of convective driving is the Bolgiano–Obukhov length,
bo vT
54 34
=
ℓ
, where
T
is the average rate of thermal energy dissipation. This length scale separates
convectively-driven scales of the flow
bo
>
ℓ
ℓ
from the range of scales where the temperature fluctuations
behave as a passive scalar
bo
<
ℓ
ℓ
. In table 1this length scale is averaged over the simulation time, normalized to
the height of the simulation volume, and recorded as
bo
ℓ
¯
. The table also includes the mean Alfvén ratio,
r
A
EE
vb
=á ñ
, the time-average of the kinetic energy divided by the magnetic energy B2
2
Eb=. In the present
numerical experiments, Navier–Stokes turbulence displays the weakest form of spatial coherence while
Boussinesq magneto-convection exhibits anisotropy with regard to the direction of gravity as well as the
occurrence of large-scale spatially-coherent structures. Additionally, a dynamical anisotropy arises because of
the presence of magnetic fields.
The positions of Lagrangian tracer particles are initialized in a homogeneous random distribution at a time
when the turbulent flow is in a statistically stationary steady state. The total number of particles in the
simulation,
n
p, is listed in table 1. We use at least a million particles for a 512
3
grid. This is a standard spatial
density of tracer particles used to describe homogeneous turbulence (e.g. [53–55]). The Lagrangian statistics we
produce have been tested and found to be well-resolved in space and time; we reproduce these statistics with half
the particles. At each time step the particle velocities are interpolated from the instantaneous Eulerian velocity
field using either a trilinear (for simulations HC and MC)or tricubic (for simulation NST)polynomial
interpolation scheme. Particle positions are calculated by numerical integration of the equations of motion
using a predictor-corrector method. For the convex hull calculations, the Lagrangian particle data is resampled
at a rate of approximately
10t
h
for simulations NST and HC. The rate of sampling for simulation MC was
smaller by a factor of 10, and this was not found to impact the dispersive results examined here. Each simulation
is run for a sufficient time that Lagrangian particle pairs have separated, on average, at least by the length of the
simulation volume. We call this time the Lagrangian crossing time,
L
C
T
, and it is listed in the table in units of the
5
New J. Phys. 19 (2017)065006 J Pratt et al
Kolmogorov time scale. Lagrangian particle pair dispersion statistics exhibit a diffusive trend near this time since
the velocity fluctuations over this time and distance exhibit low correlation.
3. Pair dispersion of Lagrangian tracer particles during homogeneous Boussinesq
convection
This section presents results for particle-pair dispersion during homogeneous Boussinesq convection for
comparison with many-particle dispersion calculated from a convex hull analysis. For an introduction to the
rich field of Lagrangian particle-pair dispersion, we refer the reader to the review of Salazar and Collins [56].In
addition to this review, several more recent works [57–59]propose new dispersion phenomenologies based on
locally ballistic dynamics, an alternative to the classical idea of turbulent diffusion exhibiting scale-dependent
diffusivity [19]. Here we briefly recall the basic argument for scaling regimes of pair dispersion. For times short
compared with the autocorrelation time of the Lagrangian velocities, the relative velocity of the particles is
approximately constant. The mean-squared separation of a pair of Lagrangian particles is therefore expected to
grow quadratically with time for a short time. This is called the ballistic or Batchelor regime. The extent of the
ballistic regime is known to depend on the initial separation of the particle pair,
0
D
, due to a finite correlation of
0
D
and the rms velocity fluctuations on this scale
v
0
D
. Recent theoretical [57–59]and experimental [60,61]
works make use of a key time scale linked to
0
D
, the initial nonlinear turnover time v
00
0
tºD D. In the inertial
range of Navier–Stokes turbulence this initial turnover time can be estimated as
v2
02v
0
t~
D
()
. For times
much larger than the autocorrelation time of the Lagrangian velocity, the velocities of a pair of Lagrangian
particles are statistically independent. The mean-squared separation of a pair of Lagrangian particles is expected
to grow linearly with time. This is typically called the diffusive regime. In between the ballistic regime and the
diffusive regime is a period of time where the mean-squared separation of particle pairs can grow cubically with
time. This is typically called the Richardson–Obukhov regime. The temporal separation of the ballistic and
Richardson–Obukhov regimes [59]can be estimated by
0
t
. Achieving a clear Richardson–Obukhov regime in
direct numerical simulations depends on the initial separation of particles as well as the size of the inertial range,
and is the subject of current ongoing research for Navier–Stokes turbulence. For this reason, and due to the
limited extent of the inertial scaling range that is expected for the Reynolds numbers we obtain, we make no
claims of observing a Richardson–Obukhov regime in the present convection simulations. We compute the
initial turnover time via a one-dimensional Eulerian kinetic energy spectrum as
kk
00
301
2
Ev
t=
-
(())
with
k2
00
p=D
. Although themoderate Reynolds numbers of the presentsimulations are far away fromvalues
where a true inertial range, devoid of influences from largest or smallest scales of the flows could be realized, for
descriptive convenience we will apply this term to the interval of time-scales between the ballistic and the
diffusive regime. Figure 2illustrates the Lagrangian particle-pair dispersion for simulations HC and MC, both
driven with homogeneous Boussinesq convection characterized by a large Bolgiano–Obukhov length. In this
figure, thin solid lines indicate Batchelor scaling ∼t
2
and diffusive scaling
t~
, around the shortest and longest
timescales, respectively. For both cases, HC and MC, we have 0h
D
, and thus 0
tt
h
. Indeed, both curves
deviate from t
2
-scaling after
0
t
. For intermediate times
tt
1
0 100
000
tt-()
, they display a phase of fast
separation which eventually levels off toward diffusive dispersion. The onset of fastpair separation in convection
Figure 2. Mean-square of the separation in the direction ofgravity for pairs ofLagrangian tracer particles, dispersing in the
hydrodynamic convection simulation HC (a)and the MHD convection simulation MC (b). Particle pairs are initially separated in the
direction of gravity by 0kol
hD= (HC)and 1.4
0kol
hD= (MC). Thin solid lines indicate: Batchelor scaling ∼t
2
(short timescales)and
diffusive scaling
t~
(long timescales). Time and length are given in units of the initial turnover time
0
t
and the Kolmogorov
microscale kol
h
, respectively.
6
New J. Phys. 19 (2017)065006 J Pratt et al
at approximately
1
0
0
t
is delayed compared to Navier–Stokes turbulence where it has been observed [59]to begin
at
tt
00
t-
(
)
. In a simulation of convection an anisotropy exists between the direction of the mean
temperature gradient and the direction perpendicular. The separation of particle pairs evolves differently in
these two directions; the separation of particle pairs can also evolve differently depending on whether the pair of
particles are initially separated in the direction of the mean temperature gradient or perpendicular to it.
During Boussinesq convection with large Bolgiano–Obukhov length the Batchelor regime for pair
separations looks similar to randomly forced hydrodynamic turbulence driven at the large scales, as shown by
e.g. Sawford [62], Yeung and Borgas [63]. During the diffusive regime, large-scale flow structures associated with
large Bolgiano–Obukhov length Boussinesq convection clearly affect the pair dispersion curve. Thedispersion
curve does not look as smooth as the result obtained from randomly forced hydrodynamic turbulence driven at
the large scales. This is not surprising because the separation of the particle pairs has reached sizes comparable to
the large-scale convective plumes. We note that although our convection simulations use quasi-periodic
boundary conditions, figure 2is not qualitatively different from figure 2 of Schumacher [41], which presents
Lagrangian dispersion during Rayleigh–Bénard convection. For pair dispersion in simulations HC and MC,
extensive averaging over different flow realizations would be necessary toachieve a perfectly smoothand
universal result, free from the influence of intermittent plumes or large-scale magnetic structures.
4. Convex hull analysis of dispersing tracer particles
4.1. Description of the convex hull calculations
We seed a number of tracer particles in the simulation volume, which produces a fixed density of tracer particles.
In simulations NST, HC, and MC the number of tracer particles and their density is based on the number needed
to produce well-resolved Lagrangian pair dispersion statistics. A convex hull analysis could potentially make use
of a significantly higher density of tracer particles. To calculate a convex hull, we select and mark a group of
Lagrangian tracer particles initially contained in a small cubic sub-volume of our simulation. The initial length
scale of the group of particles
hul
l
ℓ
is calculated as the side-length of an initial cubic sub-volume; in the limit
where the group consists of only two particles,
hul
l
ℓ
would be equivalent to
0
D
, the initial separation of a particle
pair. For the density of tracer particles in simulations NST, HC, and MC,
hul
l
ℓ
varies between 20 and 30 kol
h
. The
dependence of convex hull statistics on the initial length scale and density of the particle group is examined in
appendix.
Selection of each particle group based on the initial position of the tracer particles yields groups that contain
nearly the same number of particles, with random variation of approximately 20% based on the homogeneous
random initialization of the Lagrangian tracer particles. The average number of particles in a group,
n
pch
, listed
for simulations NST, HC, and MC in table 1, is between 24 and 214. We follow the
N
hull convex hulls (see table 1)
of the marked particle groups for the span of the simulation. The required calculation of the hulls for each time-
step is performed using the standard QuickHull algorithm [64,65], implemented in the function convhulln in
the package geometry publicly available for R, from the R Project for Statistical Computing [66,67]. The surface
area and volume of the convex hulls are obtained based on a Delaunay triangulation of the hull vertices. We stop
tracking the convex hull of a group of particles when the Lagrangian crossing time,
L
C
T
, is reached to avoid the
possibility of numerical artifacts due to the periodicity of the simulation volume.
The initial positions of particle groups could be chosen in regions of special interest in the flow, but in this
work we restrict ourselves to a homogeneous initial distribution of the groups. For each simulation the ensemble
of particle groups is initially selected to fill completely a horizontal slab. The total number of groups of particles
that we analyze using convex hulls is listed as
N
hulls
in table 1. This large number of convex hulls is more than are
required for statistical convergence of average quantities, but allows us to capture some statistically rare flow
features.
As any pair of particles separates in a turbulent flow, the particles move with the small-scale fluctuations of
the velocity field. The distance between the two particles increases monotonically in time on average, but any
specific pair of particles will produce an erratic, noisy signal. If a convex hull is defined by a very small group of
particles, then most of the particles define the surface of the convex hull. These particles on the surface of the
convex hull are called vertices of the convex hull. In the situation where most of the particles are vertices, the
convex hull, like the particle-pair distance, shrinks or grows erratically as its component tracer particles move in
the turbulent flow. The limit where groups contain only small numbers of particles is of little physical interest for
convex hull analysis, because particle pairs or particle tetrahedra already provide useful dispersion information.
In simulations NST, HC, and MC we examine the relative dynamics of larger groups of particles. If a particle
that is a vertex of the convex hull moves inward toward the center of the larger group of particles, it is unlikely
that it will remain a vertex because of the requirement of convexity. It can become an interior particle of the
convex hull. Other particles may continue to move away from the group, and the convex hull will typically
7
New J. Phys. 19 (2017)065006 J Pratt et al
continue to expand smoothly. The particles that constitute the group of vertices of the convex hull can be
exchanged frequently. This is a distinctive concept for the convex hull because it provides a contrast with more
common Lagrangian diagnostics such as particle pairs or particle tetrahedra. For statistics constructed from
particle pairs or particle tetrahedra, the same particles define the size at each point in time.
The convex hull also intrinsically links a macroscopic length scale, the size of the convex hull, with the
position of the convex hull’s geometrical center. Over this length scale, the convex hull filters out tracer particles
which disperse slower than its vertices, selecting the most efficiently dispersing members of the particle group.
4.2. Convex hull description of a group of tracer particles
A convex hull is defined by its vertices; these are the particles that dispersed the fastest in a given group of
particles. Potentially this could decouple the convex hull from the enclosed particles in two ways. The number of
vertices of the convex hull could become extremely small, or the majority of interior particles could detach from
the convex hull vertices and clump somewhere in a subregion inside the hull. In this section we devise simple
basic checks for these two scenarios.
If the particles contained in the convex hull do not spread throughout the space inside of the convex hull
evenly as it grows, the convex hull will fail to characterize the full group of particles. We use the average
difference between the geometric center of the convex hull, cvtx, and the virtual center of mass of the interior
particles contained in the convex hull,
c
in
t
, as an indicator of decoupling. This difference will not be zero,
because the particles that make up the convex hull will never fill the space perfectly evenly. Since this difference
will grow in time as the particles disperse, we compare it to a maximal extent of the convex hull at any point in
time, defined by dddd
xyz
222
12
=++()
where d
x
is the extent of the convex hull projected on the x-direction,
and d
y
and d
z
are defined similarly. Figure 3(a)shows that the average difference between the centers normalized
by the convex hull’s maximal extent,
cccd
vtx int
d
=á - ñ∣∣
. This normalized average difference between centers
does not become larger than 40% during an initial phase (0.2 LCT
t0.4
LCT)and during the subsequent
phase converges toward a quasi-constant level ranging between 15% and 20%, which is less than the standard
deviation of the coordinates of the group of tracer particles for each simulation. In figure 3(a), time is given in
terms of the Lagrangian crossing time, LCT.
The differences in the initial separation of the particle pairs in figure 2(0kol
h
D
)and the mean initial
length scale of the convex hulls (20 30
hull kol
h-
ℓ
)generate dispersion curves that reflect different ranges of
temporal and spatial scales of the underlying turbulence. Because the observable dispersion regimes and their
duration can change as a consequence of different
0
D
or
hul
l
ℓ
, a direct comparison of both figures is difficult. In
Navier–Stokes turbulence, the initial turnover time,
0
t
, has been shown [57,59]to signal the transition from the
ballistic to the inertial range of dispersion, and thus to provide a reference scale of dispersion. We therefore use
the initial turnover time
0
t
to normalize the dispersion of particles contained in a convex hull, with initial length
scale
hul
l
ℓ
. It is however not expected that this normalization can eliminate the physical differences between
isotropic Navier–Stokes and anisotropic convective systems. Moreover, it is not clear whether the universality of
0
t
extends beyond the transition from ballistic to Richardson-like dispersion.
Close examination of figures 2and 3shows that the initial phase, during which the average difference in the
centers of convex hull and interior particles increases to a maximal value, extends into the fast separation regime
of particle pair dispersion. The subsequent phase of decreasing
cd
corresponds to separation scales near to and in
Figure 3. (a)The average distance between the geometric centers of the convex hulls,
c
vtx
, and the virtual center of mass of their
interior particles, cint, divided by the convex hull size, dddd
xyz
222
1
2
=++()
.(b)Data as shown in panel (a), shifted vertically to
common initial value. Solid line: NST, dotted line: HC, dashedline: MC. Averaging is performed over convex hulls calculated for each
group of Lagrangian tracer particles and at each time. Time is given in units of (a)the Lagrangian crossing time, LCT, and (b)the initial
turnover time,
0
t
.
8
New J. Phys. 19 (2017)065006 J Pratt et al
the diffusive regime. These signatures, as well as the sharp transients evident between phases of the evolution of
cd
in figure 3, indicate a potential utility of the convex hull for studies of the turbulent inertial range.
Figure 4reveals the distribution of the group of particles within the convex hull in the z-direction. Here the
z-direction has been selected because it is the direction of the gravitational anisotropy in the convective cases;
however for one-dimensional cuts in directions other than the z-direction, similar curves result. The ratio
plotted in figure 4is the standard deviation of the particle positions,
pz,
s
, divided by the extent of the hull in the
z-direction. This ratio would be small if many of the tracer particles were to form a clump rather than spreading
throughout the interior of the convex hull. For each of the three simulations we study, however, this quantity
quickly comes to a plateau. After 0.1–0.2 LCT, i.e. the scales probed by the particles as they approach the inertial
range, the ratio no longer decreases substantially.
In all simulations the average number of vertices of the convex hulls decreases only mildly with time; this
decrease is on the order of 10% before the Lagrangian crossing time is reached. After the short initial phase up to
0
t
, the decrease in the number of convex hull vertices happens very gradually.
We conclude that on average in simulations NST, HC, and MC, the convex hulls and their interior particles
do not detach from each other in a way that would render the concept of the convex hull inappropriate for
characterizing a pre-selected group of many Lagrangian particles. Based on the measurements presented, a clear
distinction can be made between the diffusive regime and the inertial range. The correlated inertial-range
velocity fluctuations lead to changes in the relationship of convex hull vertices and interior particles. This trend
is reversed as soon as the diffusive regime is reached, largely neutralizing the differences between interior
particles and the convex hull vertices on inertial scales. This susceptibility of the convex hull to the different
characteristic regimes of ballistic, inertial-range, and diffusive turbulent transport render this diagnostic
attractive for future Lagrangian investigations of turbulence.
Apart from the ability of the convex hull to indicate different regimes of turbulent transport, the tests above
also yield information about the dynamics of the turbulent velocity field. The average displacement shown in
figure 3quantifies anisotropic differences between the dynamics of the most efficiently dispersing convex hull
vertices and the slower dispersing interior particles. On the spatial scale set by the convex hull, an anisotropic
difference of the velocity fluctuations responsible for vertex and interior tracer transport is observable as a
relative displacement of the centers of the group of interior particles and of those that define the convex hull.
Figure 3(b)which is a different representation of the data shown in figure 3(a)demonstrates this point. All three
systems have slightly different initial Lagrangian tracer configurations and, consequently, the corresponding
initial values of
cd
differ by up to 8% (MC)while NS and HC have an initial difference of approximately 1%.
Shifting the
cd
-curves of all three systems to a common initial level allows a qualitative comparison although this
simple approach can not eliminate all dynamical differences caused by varying initial tracer separations.
The increase observed for
cd
is driven by the particles that are part of the surface of the convex hull, since they
determine the geometric center of the convex hull. The relative motion of particles contained in the interior of
the hull is driven by velocity fluctuations on scales smaller than the convex hull size. Particles at opposite
locations on the surface of the convex hull will experience velocity differences on the scale of the convex hull and
therefore tend to move more rapidly apart from each other than particles in the interior of the hull, which in turn
determine the center of mass of the convex hull. Thus on time scales of
tt
00
t-
(
)⪅
a significant displacement
between the geometric center and the center of mass of a group of particles can occur, evidenced in the rapid
growth of
cd
. The relative displacement of the geometric center and the center of mass continues to grow at a
slower rate for
tt
00
t->
(
)
. This can be attributed to a finite time correlation of the velocity fluctuations on the
Figure 4. The standard deviation pz,
s
of the zcoordinates of interior particles of a convex hull divided by its extension along the
z-direction d
z
. Averaging is performed over all convex hulls in each simulation. Time is given in units of the Lagrangian crossing time,
LCT.
9
New J. Phys. 19 (2017)065006 J Pratt et al
scale of the convex hull. In addition, as time evolves and the hull grows in size, particles in the interior of the
convex hull will also experience increasing velocity fluctuations and thus some interior particles may become
particles on the surface of the hull and—vice versa—particles on the surface of the convex hull can move into its
interior due to engulfment by other particles. This process eventually leads to a decrease in the relative
displacement of the geometric center and the center of mass as the diffusive regime is approached. A noticeble
difference between the NST configuration and the convective systems HC and MC is the presence of a plateau
for the NST case between
3
0
t
and
1
6
0
t
, while for HC and MC,
cd
continues to grow during this time. The
different behavior may be caused by anisotropy in the convective flows HC and MC sustaining longer
correlations in time for velocity fluctuations in preferential directions, which does not occur for the statistically
isotropic Navier–Stokes case.
The quantity shown in figure 4measures the diffusive character of the motion of the interior particles, rather
than dynamical anisotropy. This measure exhibits a rapid transient around
0
t
from initial levels towards a first
roughly constant plateau throughout the inertial range that finally approaches the asymptotic diffusive value
around tt 100
00
t-
(
). Here, the inherently hydrodynamic simulations NST and HC display less variation
throughout the inertial range than system MC which exhibits additional flow structuring due to the presence of
magnetic field fluctuations. This brief interpretation allows for extensions, for example focusing on vertex
dynamics or a detailed direction-specific analysis by introducing spatial projections of the hulls to narrow down
the structure of the underlying anisotropic fluctuations. This will be subject of future work.
4.3. Multi-particle dispersion using convex hull analysis
Because ballistic and diffusive ranges for particle pair dispersion are typically discussed in terms of length
squared, we employ analogous measures for a group of particles and convex hulls. This is intended to make
comparison with dispersion curves as simple and direct as possible. We calculate a maximal ray rinternal to a
convex hull defined by a group of particles G:
rxxyyzzmax . 5
ij G ij ij ij
,
222
=-+-+-
Î()()() ()
By definition, the particles
i
j,that contribute to the maximum in this definition are always vertices of the convex
hull. If the group of particles densely filled a sphere, the convex hull would be the surface of the sphere, and the
maximal ray would be the diameter of the sphere. For this reason the maximal ray is sometimes also called the
diameter of a convex hull. However in this work we examine anisotropic systems where the convex hull of a
group of particles is not typically close to spherical; we opt for the more accurate former term. The susceptibility
of the maximal ray’s orientation to deformations of the convex hull can be parameterized by the RMS value of
the vertex distance from the hull’s geometrical center normalized by the average distance to the center (averages
taken over the hull vertices),Qr
rvtx
vtx
s=.IfQ0», i.e. the convex hull is close to spherical, the maximal ray
can change its direction by an arbitrary amountand much faster than the autocorrelation times ofthe underlying
turbulent fluctuations would suggest. In this case, small fluctuations of the hull radius, which can occur due to
uncorrelated small-scale fluctuations, will lead to rapid changes of orientation of the maximal ray. The maximal
ray is highly susceptible to anisotropic deformations of the convex hull. In contrast, a significant anisotropic
deformation of the convex hull, Q1, acts like a threshold for directional variation of the maximal ray and
stabilizes its orientation. This subsection focuses onquantities specific for the convex hull and their relation to
classical Lagrangian mean-square pair-separation,
02
á
D-D ñ()
.
We average the square of the maximal ray over all groups of particles; the results for each simulation are
shown in figure 5(a). This figure demonstrates that the maximal ray, although it is not tied to the same particle
pair in each tracer group, asymptotically converges to a ballistic regime signature ∼t
2
up to approximately
0
t
,
and an asymptotic diffusive regime
t~
at long times, for all systems considered. The data shown for MC does not
attain the same temporal resolution as for systems NST and MC due to a larger time step but penetrates further
into the diffusive regime.
Additional length-scale estimates can be obtained from taking appropriate powers of the normalized surface
area,
r
S4
S12
p=(())
, and the volume,
r
V34
V13
p=(())
, of the convex hulls. Averaging the length scales
produced by many different particle groups reveals dispersive behaviors that also tend to obey the ballistic and
diffusive scaling laws. A comparison of dispersion curves produced from the surface area and volume of
simulation NST are shown in figure 5(b). Similar to Lagrangian pair dispersion, the expected asymptotic scaling
laws for ballistic and diffusive regimes are approached by the surface- and volume-based distance
approximation. However, they hold over a shorter period of time than thoseshown in figure 2. Although
figure 5(b)shows dispersion curves only for simulation NST, similar results are foundfor simulations HC and
MC. A Richardson–Obukhov-like regime is not observed. Because achieving a clear Richardson–Obukhov
regime in direct numerical simulations depends on the initial separation of particles as well as the size
of the inertial range, a Richardson–Obukhov regime is not expected in our simulations. Particle filtering, an
inherent property of the selection criterion of convex hull vertices, may also contribute to the lack of a clear
10
New J. Phys. 19 (2017)065006 J Pratt et al
Richardson–Obukhov regime resulting from convex hull analysis of dispersion. During early dispersion, the
vertices of a convex hull tend to be particles that move away from the center of the hull most rapidly in the
direction radially outward from the center of the particle group; this may explain the quasi-ballistic signature
before approximately
1
6
0
t
(see figure 3(b)). As noted by Bianchi et al [50], although there is a conceptual
connection between many-particle groups and particle pairs, many-particle groups provide different
information when measuring dispersion scalings.
There is a fundamental difference between the maximal ray and the surface- or volume-based length
approximations that becomes particularly important with regard to deformations of the convex hull: the
maximal ray by definition runs along the direction of maximum extent of the convex hull. In contrast, the other
two quantities yield averaged and isotropized approximations of the length scale probed by the hull, i.e. the
radius of a reference sphere of same surface or volume. Spherical geometry is a natural first-order approximation
of a convex hull or, more precisely, the convex polyhedron, which we use as its numerical representation, since
convexity implies that the hull has no corner pointing inwards. This constraint severely restricts the complexity
of the hull’s surface structure, since any such corner vertex would turn into an interior point enclosed by the hull.
This results in an object which can mainly be deformed by flattening of the inscribed spheroid along some
direction perpendicular to the maximal ray. The convex hull is not material and therefore is not constrained by
volume conservation in incompressible flow. Although the possible length definitions do not show large
qualitative differences compared to the maximal ray, their behavior relative to each other reflects the different
responses of hull area and volume to deformations of the convex hull. This will be exploited in section 5.
5. Results: anisotropic dynamics of convex hull vertices
The relationship between the surface area, S, and volume, V, of a convex hull reveals the anisotropy of vertex
transport in a turbulent flow, which is of particular interest during convection, i.e. in the presence of coherent
velocity structures. We introduce the non-dimensional ratio
SV
23
as a direct way to quantify anisotropy.
Because a sphere minimizes the amount of surface area for a given volume, an absolute lower bound of
4
43 4.8
23
pp»()
exists for this non-dimensional surface–volume ratio. An anisotropic convex hull, e.g. a
cigar-shaped or a pancake-shaped hull, will have a higher surface to volume ratio, so the ratio gives an
impression of how non-spherical the current state of the hull is. The ratio cannot differentiate between prolate
(cigar-shaped)and oblate (pancake-shaped)convex hulls, because it approaches infinity in the limit both of zero
pancake thickness and infinite cigar length. Higher values indicate a basiclevel of anisotropic deformation.
Figure 6(a)shows the time evolution of the surface–volume ratio, averaged over all convex hulls in each
simulation. Because the particle groups consist of small numbers of particles which are randomly distributed,
they are not initially perfectly isotropic and do not evenly fill the cubic initial volumes; the resulting convex hulls
do not form either perfect cubes or perfect spheres. Thus the surface–volume ratio initially exhibits an average
value of approximately 5.6, a low value that lies between the values for perfectly spherical and perfectly cubical
volumes. In all simulations, the surface–volume ratio begins to increase around
t
t=
h
, indicating that the
convex hulls typically become stretched, i.e. anisotropic, as their particles start to disperse due to turbulent
fluctuations. In the Navier–Stokes case (NST)no global anisotropy exists in the flow. As expected, the average
surface–volume ratio remains relatively low throughout the simulation reaching its maximal value around 10 th.
Figure 5. (a)Evolution of the mean-square maximal ray rof the convex hulls in all three systems. (b)Evolution of mean-square
maximal ray X=r(solid curve), of the length basedon the hull’s surface area,
X
S4
12
p=()
(dot-dash), and of the length based on
the hull’s volume,
X
V34
13
p=()
(dash-3dot), for simulation NST, thin solid lines as in figure 2. Brackets indicate averaging over all
groups of tracer particles in a horizontal slab in each simulation volume. The symbols r
0
,S
0
, and V
0
denote the respective initial values.
11
New J. Phys. 19 (2017)065006 J Pratt et al
At long times, the average surface–volume ratio returns to approximately its initial value as uncorrelated particle
motion begins to eliminate anisotropic deformations of the convex hull. The changes in the surface–volume
ratio also slow and it approaches a flat regime related to the diffusive trend observed in figure 5at long times.
In the case of hydrodynamic Boussinesq convection (HC), the mean temperature gradient introduces a
preferential direction. We would thus straightforwardly expect higher stretching of the hulls in this direction.
However, figure 6shows that this does not take place for the convex hulls we followed; for times greater than th
only a slight increase occurs followed by a plateau phase up to
3
0t
h
. Subsequently, the average surface–volume
ratio quickly decreases below its initial value. The scale of the convective plumes in simulation HC are large and
diffuse, as reflected by the large Bolgiano–Obukhov length bo
ℓ
, the smallest scale on which the cascade of
thermal fluctuations is driven by buoyancy [68]. This large Bolgiano-Obukhov length indicates that smaller-
scale turbulent dynamics are not driven by the anisotropic influence of buoyancy. The convex hulls do not tend
to become strongly anisotropic, because the length scale of the anisotropic convection differs considerably from
the scale of the convex hulls examined. The Reynolds number of HC which is less by approximately 50% as
compared to the value of the NST system also explains why the HC simulation exhibits the lowest level of
convex-hull anisotropy.
A different behavior is observed for the magnetohydrodynamic convection (MC)simulation since larger-
scale magnetic fluctuations have a strong impact on small-scale dynamics (in contrast to a large-scale velocity
there exists no frame of reference which eliminates the magnetic field); consequently, far higher surface–volume
ratios are attained than in the other two cases. In this simulation the large-scale magnetic field fluctuations result
in strong local anisotropy of the small-scale velocity fluctuations [69–76]; the consequence is considerable
stretching of the convex hulls.
The mean alone does not characterize the full information that the convex hull analysis can provide about
anisotropy in each simulation. The shape of the probability distribution of the surface–volume ratio yields a
more comprehensive picture. If all convex hulls in a simulation were perfect spheres, the distribution of the
surface–volume ratio would be a delta function. However, the distributions show a strong dependence on the
type of turbulence indicated by the values of distribution mean, μ, and standard deviation, σ, given in the
caption of figure 6. In the hydrodynamic convection case the distribution is the narrowest with the lowest mean
indicating the highest level of anisotropy, followed by NST and MC. The significant hull anisotropy observed for
system MC is a clear indication of the additional anisotropy imposed by the slowly evolving large-scale magnetic
field fluctuations on the smaller-scale velocity fluctuations. These results are not surprising and consistent with
the data given in figure 6(a). In addition, figure 6(b)shows the centered and normalized distributions of the
surface–volume ratio for each of our three simulations after each set of convex hulls has evolved for 20 th. All
distributions collapse on a positively skewed functional shape suggesting a general characteristic of convex hull
deformation common to all three turbulent systems.
The surface–volume ratio varies spatially in each simulation. The time evolution of this ratio for a single
convex hull in simulation MC is illustrated in figure 7. At early times, the surface–volume ratio for this individual
hull grows to considerably exceed the mean, indicating that this hull is more stretched than the average convex
hull of this ensemble. This surface–volume ratio also exhibits rapid changes in time. For example, during the
Figure 6. The time evolution of the convex hull’s surface area, S, divided by the 2/3 root of the volume, V.In(a)the evolution of this
non-dimensional surface–volume ratio is averaged over all convex hulls in each simulation. In (b)the probability distribution
function, P,of
SV
23
ms-()
is shown at time 20
t
h
,μdenoting the mean and σthe standard deviation of the respective
distribution. The tuples ,ms(
)
are NST:(6.8, 0.7), HC:(6.0, 0.3), MC:(8.5, 1.6).
12
New J. Phys. 19 (2017)065006 J Pratt et al
period between approximately 5 thand 10 ththis hull goes from a more anisotropic form than average to a
considerably less anisotropic form.
In figure 7, the surface–volume ratio is also shown as a contour plot for the set of convex hulls that fill a
horizontal slab of simulation MC. Dark areas represent regions where convex hulls have grown with significant
anisotropy. High spatial intermittency is also noticeable, with areas of large anisotropy bordering areas that grow
more isotropically. This pattern of anisotropy remains similar over a long period of time, reflecting the strong
influence of the initial configuration of the flow on local dispersion. Although we examine a small number of
simulations, the non-dimensional surface–volume ratio that we introduce is clearly capable of revealing aspects
of local anisotropy in turbulent flows.
6. Results: extreme-value statistics of turbulent particle dispersion
The vertices of a convex hull are the particles that disperse fastest among a given group of particles, and the
maximal ray defines a maximal dispersion of all particle pairs within the group. Thus the use of the convex hull
evokes concepts from extreme value theory [77,78]. The most widely encountered distribution in extreme value
theory, the Gumbel distribution [79,80], has been frequently employed for climate modeling, including extreme
rainfall and flooding [81–85], extreme winds [86], avalanches [87], and earthquakes [88]. The Gumbel
distribution has also been found to reasonably characterize the density fluctuations within galaxies [89–91]and
in certain areas of tokamaks [92–94], binding energies in liquids [95], as well as turbulent fluctuations [96,97].
The cumulative distribution function Ffor the Gumbel case has the well-known form:
Fx xexp exp , 6mb=---() ( ( ( ) )) ()
where the location parameter μgives the mode of the distribution, βis commonly called the scale parameter,
and the median of the distribution is ln ln 2
m
b-(())
. Because extreme value theory typically develops as an
asymptotic theory for sample sizes
n
~¥
, convex hulls with large numbers nof particles facilitate the
exploitation of extreme value theory results.
We examine the square-length of the maximal ray with extreme value theory, and this choice is crucial. The
square-length of the maximal ray is a fundamental scalar commonly associated with dispersion, and thus the
most natural physical quantity to consider. The square-length of the maximal ray is also consistent with a simple
model of Gaussian displacements. No rigid upper limit exists for the square-length of the maximal ray, and thus
the Gumbel distribution is the case that would be anticipated from extreme value theory.
Because Lagrangian tracer particles move in a flow with a finite correlation in space and time, their motions
are not independent. The number of particles in each group is also limited in these numerical experiments.
Despite these limitations, we find that the shape of the cumulative distribution function of the square of the
maximal ray is suggestive of a Gumbel distribution. This observation holds at each point in time, regardless of
whether the particle groups sampled are in the ballistic regime, diffusive regime, or a transitional period of
dispersion. A Gumbel distribution describes the results well, regardless of the initial length scale of the convex
hull, and the initial density of particles, for the range
4
64
hull
kol kol
hh<<ℓthat we have tested (see appendix).
Figure 7. (Left)a comparison of the non-dimensional surface–volume ratio between the convex hull of asingle arbitrarily chosen
group of tracer particles and the average, in the simulation MC. (Right)a contour plot that shows a horizontal slab filled with convex
hulls in simulation MC, at a late time in the simulation. Darker colors represent higher values of the surface–volume ratio. The colors
are shown at the initial positions of the convex hulls, and each pixel approximately represents the initial volume of a convex hull.
13
New J. Phys. 19 (2017)065006 J Pratt et al
This suggests that the Gumbel distribution might provide an effective description of the probability of extremes
of turbulent dispersion. The location and scale parameters can be different for different
hul
l
ℓ
, and at different
times in the dispersion process, although a Gumbel distribution is recovered at each time.
In addition, we consider a cumulative distribution function constructed from data at all times throughout
the evolution of the convex hulls, as shown in figure 8. Using data from all times is a reasonable choice that
produces a single form of the cumulative distribution function relevant to the entire simulation. From the
perspective of the simple model of Gaussian displacement, noted above, that pragmatic choice actualizes a
distribution of values of the scale parameter. Such a possibility is well known in related, but physically distinct,
studies of turbulence [98]. Figure 8(a)shows that the distribution of square-length of the maximal ray is fit well
with a Gumbel distribution when physically distinct directions, perpendicular and parallel to gravity, are
considered individually in the magnetoconvection simulation MC. We found in section 5that the convex hulls
in simulation MC become highly anisotropic on average. Thus the fact that a Gumbel distribution with different
location and scale parameters accurately describes the extremes of dispersion in both of these physically distinct
directions is a new and significant physical observation.
In figure 8(a), we observe an ordering between the scale parameter obtained for the direction perpendicular
to gravity and the direction parallel to gravity; the value of the scale parameter is larger in the direction parallel to
gravity. Figure 8(b)compares the cumulative distribution functions of the square-length of the full maximal ray
of the convex hull in each simulation, and again they demonstrate the linear behavior expected of
F
l
nln-(())
for
the Gumbel distribution. When the
F
l
nln-(())
is fit using linear regression, the value of the scale parameters
are:
0.17
HC
b=
,
0.40
NST
b=
, and
0.65
MC
b=
. In section 5, ordering these simulations according to the least
anisotropic to the most anisotropic simulation produced: HC, NST, MC. We thus conjecture from the results in
figures 8(a)and (b)that faster dispersion linked to anisotropy will lead to a higher value of the scale parameter.
7. Discussion
We have shown that the convex hull can be used to characterize many-particle dispersion in turbulent flows, and
can reproduce scalings similarto particle-pair, and other multi-particle Lagrangian statistics.The convex hull
allows us to extract dispersion behaviors that produce clear scalings from groups of tracer particles that are
significantly larger than have been typically examined by multi-particle statistics. We have examined particle
dispersion using convex hulls across three types of physically distinct turbulence simulations, including
Navier–Stokes turbulence, Boussinesq convection, and MHD Boussinesq convection. In each of the simulations
that we consider, we have shown that the convex hull describes well the dynamics of the entire group of particles.
In addition, these tests yield further information about the turbulent velocity field by quantifying the dynamical
differences between interior particles and convex hull vertices. Dispersion curves produced using the maximal
ray of the convex hull, the surface area of the convex hull, and the volume of the convex hull produce ballistic and
diffusive scalings, which can be compared with particle-pair dispersion curves. Although the convex hull has
Figure 8. The log negative log of the cumulative distribution function, F, of the square of the maximal ray of the group of particles
defined in equation (5). Panel (a)shows the cumulative distribution function of the square of the maximal ray in the directions
perpendicular and parallel to gravity from the MHD convection simulation MC. In (b)shows the cumulative distribution function of
the square of the maximal ray for each simulation. Foreach cumulative distribution function shown, a line (solid black line)fits the
natural log of the negative natural log of Fwell.
14
New J. Phys. 19 (2017)065006 J Pratt et al
been used to calculate volumes occupied by particles in some specialized contexts [21,27], this is the first time
that the convex hull of the positions of Lagrangian tracer particles has been used as a fundamental diagnostic to
obtain Lagrangian statistics of multi-particle dispersion in homogeneous turbulent flows.
In addition, we have explored the convex hull’s fundamental link to extreme value statistics. We have
discussed that the convex hull provides new information about extremes of dispersion that standard
multi-particle statistics cannot. Convex hulls calculated from large numbers of particles provide an ideal
application for extreme value theory, an asymptotic theory for large samples. Predictions based on extreme value
theory are of practical use for studies of contaminants or of energetic particles, where questions about maximal
dispersion are critical. Experimentally it may be simpler to track the convex hull of a large number of particles
than to track all the particles in the group individually. We show that the distribution of the square length of the
maximal ray of the convex hull is the Gumbel case of generalized extreme value distributions. In addition we
show that for a system that is anisotropic because of MHD convection, the maximal ray in each physically
distinct direction is described well by the Gumbel distribution. Because the Gumbel distribution has been
successful in predicting avalanches, extreme rainfall, and extreme winds, thisnontrivial new observation will
provide new physical intuition for modeling anomalous dispersion.
In a second application of the convex hull analysis, we exploit the relationship between convex hull surface
area and volume to examine the degree of anisotropy present in a turbulent convective flow. Our results reveal
the extent of spatial variation of anisotropy. Moreover, this quantity also exhibits a probability distribution that
has the same universal shape for all three considered physical systems. Convex hull analysis can easily isolate
dispersive characteristics in any local region of interest, for example a region where a magnetic structure, or
strong convective plume is present. Used in this way, they provide a versatile supplement to standard Lagrangian
multi-particle statistics in complex turbulent flows. Because of these advantages, further investigation of the
convex hull to analyze many-particle turbulent dispersion is justified.
Acknowledgments
We thank Luca Biferale for his helpful comments on this work. The research leading to these results has received
funding from the European Research Council under the European Union’s Seventh Framework (FP7/2007-
2013)/ERC grant agreement no. 320478. This work has also been supported by the Max-Planck Society in the
framework of the Inter-institutional Research Initiative ‘Turbulent Transport and Ion Heating, Reconnection
and Electron Acceleration in Solar and Fusion Plasmas’of the MPI for Solar System Research, Katlenburg-
Lindau, and the Institute for Plasma Physics, Garching (project MIFIF-A-AERO8047). Original simulations
were performed on theVIP, VIZ, and HYDRA computersystems at the Rechenzentrum Garching of theMax
Planck Society. Additional calculations were performed on the Konrad and Gottfried computer systems of the
Norddeutsche Verbund zur Förderung des Hoch- und Höchstleistungsrechnens (HLRN). NWW acknowledges
travel support from KLIMAFORSK project number 229754 and the London Mathematical Laboratory, and
Office of Naval Research NICOP grant NICOP—N62909-15-1-N143 at Warwick and Potsdam.
Appendix. Dependence of convex hull statistics on initial size and density
Because the initial separation between a pair of tracer particles affects two-particle dispersion, the initial length
scale
hul
l
ℓ
of a group of tracer particles may also play a role in a convex hull analysis of many particle dispersion.
The initial density of tracer particles clearly also is significant for dispersion, because this directly determines the
resolution of the convex hull surface area and volume. For simulations NST, HC, and MC, the initial density of
tracer particles is between 0.001 and
0
.02
3
particles
kol
h
. This density severely limits the
hul
l
ℓ
that can be explored
in these simulations. The initial size of the particle groups is chosen to be, 20 30
hull
kol kol
hhℓ, and groups
with substantially smaller
hul
l
ℓ
clearly do not contain enough particles for the convex hull to be adequately
resolved. For a point of comparison, the particle groups examined by Bianchi et al [50]are significantly smaller
and more dense; they have initial length scale of hull kol
h=
ℓ
which contains 2000 tracer particles.
In order to systematically test the convex hull analysis of dispersion for a range of initial sizes, we perform a
test simulation of forced homogeneous isotropic Navier–Stokes turbulence, similar to simulation NST, in which
tracer particles are initialized in groups with given
hul
l
ℓ
, and at two different fixed particle densities. The first
density, 0.005 3
particles
low kol
r
h=, is selected to be similar to the tracer particle density in simulations NST,
HC, and MC, in order to examine
hul
l
ℓ
both larger and smaller those examined in these simulations. At this
density we examine groups of tracer particles with eight initial sizes between
1
464
hull
kol kol
hhℓ. The second
density
11.5
3
particles
hi kol
r
h=
is significantly higher, so that we can examine groups of particles with four
smaller initial sizes between
4
14
hull
kol kol
hhℓ.
15
New J. Phys. 19 (2017)065006 J Pratt et al
Regardless of
hul
l
ℓ
and initial particle density, the trends evident for the convex hull validity diagnostic in
figures 3and 4are recovered. The time at which thenormalized average difference between centers
cccd
vtx int
d
=á - ñ∣∣
reaches a peak appears to be approximately independent of both the
hul
l
ℓ
and density. For
larger
hul
l
ℓ
the growth of
cd
begins earlier, although the time at which
cd
begins to grow is not directly relatable
to
0
t
. The decrease to a plateau, evident in figure 4, is larger for groups of particles with larger
hul
l
ℓ
,atfixed
particle density for both densities tested.
Aside from these diagnostics, we consider the growth of the difference between the geometric center of the
convex hull and the center of mass of the group of tracer particles, normalized by the initial length scale of the
convex hull
cc
hullvtx int
á
-ñℓ∣∣
. This is not useful for examining whether the convex hull describes the group of
particles well, because unlike the maximal extent d,
hul
l
ℓ
describes the initial state of the particle group and does
not change in time; however the difference in these centers provides a new quantity linked to the dispersion of
many tracer particles. The evolution of this quantity is shown in figure 9. In the figure, the magnitude of the
difference in the centers is clearly linked to
hul
l
ℓ
, with larger
hul
l
ℓ
leading to smaller values of
cc
hullvtx int
á
-ñℓ∣∣
throughout dispersion. However the shape of the evolution curves for the difference in the centers appears to be
approximately independent of
hul
l
ℓ
. During the inertial range of time scales, the difference in centers grows
linearly with time. Explaining this interesting scaling result will be the subject of future work.
The evolution of the maximal ray rintroduced in section 4.3 universally exhibits a clear diffusive regime for
all
hul
l
ℓ
and densities tested. This is illustrated by figure 10. We do not expect perfect ballistic scaling of the
maximal ray, because unlike a particle pair, the particles that determine the maximal ray of a convex hull can be
exchanged. Despite this, we do observe ascaling reminiscent of ballistic behavior forall
hul
l
ℓ
and densities tested;
this likely indicates that vertex exchange is not a dominant effect during this early regime of dispersion. The slope
of the dispersion curves during transitional regimes between ballistic and diffusive appears to be dependent on
the initial length scale
hul
l
ℓ
. This is unsurprising because it is a well-known result for particle-pairs.
Figure 9. Evolution of the difference between the geometric center of the convex hull and the center of mass of the group of tracer
particles, normalized by the initial length scale of the convex hull cc
hullvtx int
á
-ñℓ∣∣
for (a)groups of particles with four initial sizes
and particle density h
i
r
, and (b)groups of particles with eight initial sizes and particle density
low
r
. The initial size
hul
l
ℓ
is labeled in
units of kol
h
. A thin solid line with slope1 is indicated during time scales associated with theinertial range.
Figure 10. Evolution of the mean-square maximal ray rfor groups of particles with (a)four initial sizes and particle density h
i
r
, and (b)
eight initial sizes and particle density
low
r
. The initial size
hul
l
ℓ
is labeled in units of kol
h
.
16
New J. Phys. 19 (2017)065006 J Pratt et al
References
[1]Taylor G I 1922 Proc. Lond. Math. Soc. 20 196
[2]Falkovich G, Gawȩdzki K and Vergassola M 2001 Rev. Mod. Phys. 73 913
[3]Hackl J, Yeung P and Sawford B 2011 Phys. Fluids 23 065103
[4]Sawford B, Pope S and Yeung P 2013 Phys. Fluids 25 055101
[5]Lüthi B, Ott S, Berg J and Mann J 2007 J. Turbul. 8N45
[6]Toschi F and Bodenschatz E 2009 Annu. Rev. Fluid Mech. 41 375
[7]Xu H, Ouellette N T and Bodenschatz E 2008 New J. Phys. 10 013012
[8]LaCasce J 2008 Prog. Oceanogr. 77 1
[9]Yeung P K 2002 Annu. Rev. Fluid Mech. 34 115
[10]Pumir A, Shraiman B I and Chertkov M 2000 Phys. Rev. Lett. 85 5324
[11]Lin J, Brunner D, Gerbig C, Stohl A, LuharA and Webley P 2013 Lagrangian Modeling of the Atmosphere vol 200 (New York: Wiley)
[12]Pratt J, Busse A and Müller W-C 2013 Astron. Astrophys. 557 A76
[13]Mazzitelli I, Fornarelli F, Lanotte A and Oresta P 2014 Phys. Fluids 26 055110
[14]Maeder A, Georgy C and Meynet G 2008 Astron. Astrophys. 479 L37
[15]Brun A S, Miesch M S and Toomre J 2011 Astrophys. J. 742 79
[16]Leprovost N and Kim E-J 2006 Astron. Astrophys. 456 617
[17]Efron B 1965 Biometrika 52 331
[18]Gifford F Jr 1957 J. Atmos. Sci. 14 410
[19]Richardson L F 1926 Proc. R. Soc. Lond. Ser. A110 709
[20]Elder J 1959 J. Fluid Mech. 5544
[21]Stocker R and Imberger J 2003 Limnol. Oceanogr. 48 971
[22]Spencer D, Lemckert C J, Yu Y, GustafsonJ, Lee S and Zhang H 2014 J. Coast. Res. 70 29–34
[23]Schmeja S and Klessen R 2006 Astron. Astrophys. 449 151–9
[24]Áingeles Serna M, Bermudez A, Casado R and Kulakowski P 2011 Proc. 7th Int. Conf. on Intelligent Sensors, Sensor Networks and
Information Processing pp 419–24
[25]Li S et al 2008 Bioinformatics 24 i182
[26]Millett K, Rawdon E, Stasiak A and Sułkowska J 2013 Biochem. Soc. Trans. 41 533
[27]Dietzel M and Sommerfeld M 2013 Powder Technol. 250 122
[28]Cornwell W K, Schwilk D W and Ackerly D D 2006 Ecology 87 1465
[29]Randon-Furling J, Majumdar S N and Comtet A 2009 Phys. Rev. Lett. 103 140602
[30]Majumdar S, Comtet A and Randon-Furling J 2010 J. Stat. Phys. 138 955
[31]LukovićM, Geisel T and Eule S 2013 New J. Phys. 15 063034
[32]van der Zanden H B et al 2013 Mar. Ecol. Prog. Ser. 476 237
[33]Dumonteil E, Majumdar S N, Rosso A and Zoia A 2013 Proc. Natl Acad. Sci. USA 110 4239
[34]Collyer M L, Hall M E, Smith M D and Hoagstrom C W 2015 Copeia 2015 181
[35]Avellaneda M and Weinan E 1995 Commun. Math. Phys. 172 13
[36]Bertoin J 2001 Levy Processes (Berlin: Springer)pp 267–79
[37]Chupeau M, Bénichou O and Majumdar S N 2015 Phys. Rev. E91 050104
[38]Busse A, Müller W-C, Homann H and Grauer R 2007 Phys. Plasmas 14 122303
[39]Homann H, Ponty Y, Krstulovic G and Grauer R 2014 New J. Phys. 16 075014
[40]Schumacher J 2009 Phys. Rev. E79 056301
[41]Schumacher J 2008 Phys. Rev. Lett. 100 134502
[42]Forster C, Stohl A and Seibert P 2007 J. Appl. Meteorol. Climatol. 46 403
[43]Busse A 2009 PhD Thesis Universität Bayreuth
[44]Moll R, Graham J P, Pratt J, Cameron R H, Müller W-C and Schüssler M 2011 Astrophys. J. 736 36
[45]Kurihara Y 1965 Mon. Weather Rev. 93 33
[46]Williamson J H 1980 J. Comput. Phys. 35 48
[47]Gibert M, Pabiou H, Chilla F and Castaing B 2006 Phys. Rev. Lett. 96 084501
[48]Pope S B 2000 Turbulent Flows (Cambridge: Cambridge University Press)
[49]Marino R, Baerenzung J, Mininni P, Pouquet A, Rorai C, Rosenberg D and Stawarz J 2015 Direct and Large-Eddy Simulation IX (Berlin:
Springer)pp 549–59
[50]Bianchi S, Biferale L, Celani A and Cencini M 2016 Eur. J Mech. B55 324
[51]Calzavarini E, Lohse D, Toschi F and Tripiccione R 2005 Phys. Fluids 17 055107
[52]Eswaran V and Pope S 1988 Comput. Fluids 16 257
[53]Biferale L, Boffetta G, Celani A, Devenish B, Lanotte A and Toschi F 2005 Phys. Fluids 17 115101
[54]Homann H, Grauer R, Busse A and Müller W-C 2007 J. Plasma Phys. 73 821
[55]Busse A, Müller W-C, Homann H and Grauer R 2007 Phys. Plasmas 14 122303
[56]Salazar J P and Collins L R 2009 Annu. Rev. Fluid Mech. 41 405
[57]Bourgoin M 2015 J. Fluid Mech. 772 678
[58]Thalabard S, Krstulovic G and Bec J 2014 J. Fluid Mech. 755 R4
[59]Bitane R, Homann H and Bec J 2012 Phys. Rev. E86 045302(R)
[60]Ouellette N T, Xu H, Bourgoin M and Bodenschatz E 2006 New J. Phys. 8109
[61]Bourgoin M, Ouellette N T, Xu H, Berg J and Bodenschatz E 2006 Science 311 835–8
[62]Sawford B 2001 Annu. Rev. Fluid Mech. 33 289
[63]Yeung P and Borgas M S 2004 J. Fluid Mech. 503 93
[64]Barber C B, Dobkin D P and Huhdanpaa H 1996 ACM Trans. Math. Softw. 22 469
[65]Barber C B, Dobkin D P and Huhdanpaa H 2013 Qhull: Quickhull algorithm for computing the convex hull ASCL Code Record
ascl:1304.016
[66]Ihaka R and Gentleman R 1996 J. Comput. Graph. Stat. 5299
[67]R Core Team 2014 R: A Language and Environment for Statistical Computing (Vienna: R Foundation for Statistical Computing)
[68]Biskamp D 2003 Magnetohydrodynamic Turbulence (Cambridge: Cambridge University Press)
17
New J. Phys. 19 (2017)065006 J Pratt et al
[69]Grappin R and Müller W-C 2010 Phys. Rev. E82 026406
[70]Verdini A, Grappin R, Hellinger P, Landi S and Müller W C 2015Astrophys. J. 804 119
[71]Matthaeus W H, Ghosh S, Oughton S and Roberts D A 1996 J. Geophys. Res. 101 7619
[72]Cho J and Vishniac E T 2000 Astrophys. J. 539 273
[73]Chandran B D 2008 Astrophys. J. 685 646
[74]Montgomery D and Turner L 1981 Phys. Fluids 24 825
[75]Goldreich P and Sridhar S 1995 Astrophys. J. 438 763
[76]Boldyrev S 2006 Phys. Rev. Lett. 96 115002
[77]Castillo E, Hadi A S, Balakrishnan N and Sarabia J-M 2005 Extreme Value and Related Models with Applications in Engineering and
Science (Hoboken, NJ: Wiley)
[78]Majumdar S N, Comtet A and Randon-Furling J 2010 J. Stat. Phys. 138 955
[79]Bramwell S, Christensen K, Fortin J, Holdsworth P, Jensen H, Lise S, López J, Nicodemi M, Pinton J and Sellitto M 2000 Phys. Rev. Lett.
84 3744
[80]Gumbel E 1958 Statistics of Extremes (New York: Columbia University Press)
[81]Hirabayashi Y, Mahendran R, Koirala S, Konoshima L, Yamazaki D, Watanabe S, Kim H and Kanae S 2013 Nat. Clim. Change 3816
[82]Borga M, Vezzani C and Dalla Fontana G 2005 Nat. Hazards 36 221
[83]Koutsoyiannis D 2004 Hydrol. Sci. J. 49 610
[84]Coles S, Pericchi L R and Sisson S 2003 J Hydrol. 273 35
[85]Yue S 2000 Adv. Water Resour. 24 179
[86]Kang D, Ko K and Huh J 2015 Energy 86 51
[87]Schweizer J, Mitterer C and Stoffel L 2009 Cold Reg. Sci. Technol. 59 234
[88]Pisarenko V, Sornette A, Sornette D and Rodkin M 2014 Pure Appl. Geophys. 171 1599
[89]Antal T, Labini F S, Vasilyev N L and Baryshev Y V 2009 Europhys. Lett. 88 59001
[90]Waizmann J-C, Ettori S and Moscardini L 2012 Mon. Not. R. Astron. Soc. 420 1754
[91]Chongchitnan S and Silk J 2012 Phys. Rev. D85 063508
[92]Hnat B et al 2008 Nucl. Fusion 48 085009
[93]Anderson J and Kim E-J 2009 Plasma Phys. Control. Fusion 52 012001
[94]Graves J, Horacek J, Pitts R and HopcraftK 2005 Plasma Phys. Control. Fusion 47 L1
[95]Chempath S, Pratt L R and Paulaitis M E 2010 Chem. Phys. Lett. 487 24
[96]Noullez A and Pinton J-F 2002 Eur. Phys. J. B28 231
[97]Dahlstedt K and Jensen H J 2001 J. Phys. A: Math. Gen. 34 11193
[98]Castaing B, Gagne Y and Hopfinger E 1990 Physica D46 177
18
New J. Phys. 19 (2017)065006 J Pratt et al
