# Orientation of Finite Reynolds Number Anisotropic Particles Settling in
Turbulence

Anubhab Roy  Stefan Kramel  Udayshankar Menon  Greg A. Voth  Donald L. Koch

###### Abstract

We present experimental and computational results for the orientation
distributions of slender fibers and ramified particles settling in an
isotropic turbulent flow. The rotational dynamics of the particles is modeled
using a slender-body theory that includes the inertial torque due to
sedimentation that tends to rotate the particles toward a broadside
orientation. The particles are assumed to rotate due to viscous forces
associated with the turbulent velocity gradients occurring on the particle
length scale. In the simulations, the turbulence is obtained from a stochastic
model of the velocity gradient in a Lagrangian reference frame. In the
experiments, the turbulence is generated by active jets in a vertical water
tunnel. It is well known that axisymmetric particles rotate according to
Jeffery’s solution for the rotation of a spheroidal particle if one adopts an
appropriate effective aspect ratio. We show that the same result applies to a
ramified particle consisting of three coplanar fibers connected with equal
angles at a central point which rotates like a thin oblate spheroid. The
orientation statistics can be quantified with a single non-dimensional
parameter, the settling factor $S_{F}$, defined as the ratio of rotations due
to sedimentation and turbulent shear. For low values of $S_{F}$, we observe
nearly isotropically oriented particles, whereas particles become strongly
aligned near the horizontal plane for high values of $S_{F}$. The variance of
the angle away from horizontal scales as $S_{F}^{-2}$ for $S_{F}\gg 1$, but
the orientation distribution is non-Gaussian due to turbulent intermittency in
this limit.

###### keywords:

keyword one , keyword two

###### PACS:

0000 , 1111

###### MSC:

0000 , 1111

[inst1]organization=Department of Applied Mechanics,addressline=Indian
Institute of Technology Madras, city=Chennai, postcode=600036, state=Tamil
Nadu, country=India

[inst2]organization=Department of Physics,addressline=Wesleyan University,
city=Middletown, postcode=06459, state=Connecticut, country=USA

[inst3]organization=Smith School of Chemical and Biomolecular
Engineering,addressline=Cornell University, city=Ithaca, postcode=14853,
state=NY, country=USA

##  1 Introduction

Sedimentation of non-spherical particles in turbulent flows occurs in many
natural situations and has important consequences for a wide range of
engineering applications. Mixed-phase cloud systems, such as cirrus clouds,
consist of sedimenting ice crystals whose orientation distributions critically
affect global climate models. Recently lidar polarization measurements have
focussed on distinguishing droplet laden clouds from icy clouds, especially
clouds that are dominated by horizontally oriented crystals [50, 68]. A
vertically pointed Doppler lidar observes ‘mirror-like’ specular reflections
from horizontally aligned ice crystals and thus special care needs to be taken
to distinguish them from water clouds as both produce low values of the
depolarization ratio [26]. Particle shape also plays an important role in
pneumatic conveying and fluidized bed risers and current models do not account
for behavior of highly non-spherical particles such as mica flakes [24]. In
this paper we probe the competition between inertial torques due to
sedimentation that align particles and randomizing torques due to turbulent
shear in determining the orientation distribution and translational motion of
high aspect ratio particles settling in homogeneous, isotropic turbulence.
Experimental observations of ramified particles settling in a vertical
turbulent water column are complemented by theory and stochastic simulations
based on a slender-body description of the particles.

Since the orientations of small particles sedimenting in turbulence are only
affected by the nearly universal inertial and dissipation range of the flow
and not by the large scales, the dynamics of these particles is nearly the
same in many different turbulent environments. In the simplified case of
spheroids sedimenting in homogeneous, isotropic turbulence, there are five
non-dimensional parameters necessary to specify the problem. The turbulence
can be characterized by its Taylor Reynolds number,
$\textit{Re}_{\lambda}=\sqrt{15}\left(\mathcal{L}/\eta\right)^{4/3}$ where
$\mathcal{L}$ is the integral scale,
$\eta=\left(\nu^{3}/\epsilon\right)^{1/4}$ is the Kolmogorov length, $\nu$ is
the kinematic viscosity, and $\epsilon$ is the energy dissipation rate per
unit mass. A spheroid is an axisymmetric ellipsoid that can be characterized
by three non-dimensional parameters: its maximum dimension compared with the
Kolmogorov length, $L/\eta$, an aspect ratio $\kappa$ of the symmetry axis
length to a perpendicular length, and the relative density of the particle to
the fluid, $\rho_{p}/\rho_{f}$. The fifth parameter characterizes the
importance of gravity. This can be quantified by the ratio of terminal
particle velocity $W$ to the Kolmogorov velocity $u_{\eta}=(\nu\eta)^{1/4}$,
the so-called settling parameter $\textit{Sv}=W/u_{\eta}$ [4, 60, 27]. We will
see that the alignment of non-spherical particles is determined by a settling
parameter $S_{F}$ defined as the ratio of the rate of change of particle
orientation due to sedimentation induced inertia to that caused by turbulent
velocity gradients. For small fibers in the slender body limit, $S_{F}$ is
proportional to $\textit{Sv}^{2}$ and has a weak logarithmic dependence on the
particle aspect ratio.

The most accessible case of small and neutrally buoyant, non-spherical
particles in turbulence has been studied extensively in simulations [57, 69,
54] and experiments [52, 43, 37, 23] and a review is given by Voth and Soldati
[65]. The particle dynamics depend only on particle shape and possibly
Reynolds number, and preferential alignment with the local velocity gradients
results in reduced tumbling rates compared to randomly oriented particles.
This is also true for particles with lengths in the inertial range, however,
the alignment and tumbling rates of large particles depend on the coarse
grained velocity gradient at the scale of the particle [57, 53]. This
observation from experiments and direct-numerical simulations motivates our
assumption that the rotation of thin particles caused by turbulent shearing
motions can be modeled using Jeffery’s solution [29] for particle rotation
provided that a velocity gradient on the scale of the particle is used.

The larger the particle size, the more pronounced the effects of inertia when
particles and fluid are not density matched. This can alter the particle
dynamics drastically, even when gravity is neglected. Direct numerical
simulations of small, heavy ellipsoidal particles [47, 74, 42, 8] in turbulent
channel flows have shown that the preferential orientation changes non-
trivially with increasing particle inertia, especially in the near-wall
regions, which can have important consequences for particle deposition. Sabban
et al. [56] measured the translational and rotational motion of fibers in
homogeneous, isotropic turbulent air under conditions of small particle
Reynolds number and finite particle Stokes number. Most of these studies have
focused on small particles and ignored external forces.

In most physical situations, gravity cannot be ignored and heavy particles
will sediment. This again leads to very different particle dynamics and
alignment, depending on the environment, whether it is quiescent or turbulent.
Simulations [22, 72, 63] and experiments [5, 48] in quiescent fluids have
revealed great insight into the drag force and inertial torques on
sedimenting, non-spherical particles. There has been significant interest in
studying settling dynamics of ice crystals and phytoplanktons, two canonical
examples of non-spherical particles. Experimental measurements of fall
velocities and trajectories for rimed [75] and unrimed [30] ice crystals have
found complex falling trajectories depending on the Reynolds number and moment
of inertia. Ardekani et al. [2] analyzed the clustering and preferential
sampling of sedimenting phytoplanktons, modelled as inertia-less prolate
spheroids, in homogeneous isotropic turbulence. The majority of studies on
non-spherical particles in turbulence in the fluid dynamics community have
ignored the role of the inertial torque [38, 61, 8, 73, 20]. The inertial
torque determines the horizontal alignment of ice crystals settling in a
turbulent background flow, a fact that has been acknowledged in the
atmospheric science literature [9, 34, 49, 25, 51]. Inertial torques cause
slender bodies to sediment with a preferential orientation, where the long
axis is perpendicular to gravity. This is a stable orientation at low and
intermediate particle Reynolds numbers, but will eventually become unstable at
large Reynolds numbers and the particle motion becomes complex Ern et al.
[16].

The complexity of the problem becomes clear when in addition to the underlying
turbulence, the particle orientation has to be considered. Therefore, it is no
surprise that many studies focused on neutrally buoyant particles and
neglected the effects of particle inertia. In many situations however,
particles are not neutrally buoyant, there exists a density difference between
them and the fluid, and this can drastically change the dynamics and
interactions of these particles. For small particles, inertia can still be
neglected, but it becomes increasingly important with increasing particle
size. Compared to spherical particles, where inertia only affects transport
and causes enhanced sedimentation, the so-called sweeping effect ([66]), non-
spherical particles can show very different orientation distributions and
inertia can alter their preferential alignment. Direct numerical simulations
(DNS) of small, heavy ellipsoidal particles ([47, 74, 8], [42]) in turbulent
channel flows have shown that the preferential orientation changes non-
trivially with increasing particle inertia, especially in the near-wall
regions, which has important consequences for particle deposition. Most of
these studies have focused on small particles and ignored external forces.

In addition to turbulence and particle inertia, external forces such as
gravity can have a pronounced effect on particle motion. Particles with a
larger density than the fluid will sediment under the influence of gravity.
DNS of the flow field around fibers and disk-like particles ([22]) have
revealed great insight into the vast parameter space, but due to the inherent
complexity of the problem, many studies have been forced to ignore turbulence
and focus on the slightly easier task of understanding sedimentation in
quiescent fluid first. Experimental observations of free falling, non-
spherical particles have shown that they do not fall straight and do not have
random orientation distributions, but the motion depends on the particle
Reynolds number and shape [70], [75], [18] and [30]. The torque induced on
thin cylinders or prolate spheroids in this Reynolds number regime causes the
body to rotate into a stable position with its symmetry axis aligned
horizontally. This effect has been studied theoretically by [10, 39, 31], and
experimentally by [28, 5]. Only during the last decades has it been possible
to study the influence of turbulence on sedimenting non-spherical particles
([49, 48, 38]). This is of special interest to the atmospheric research
community, where prolate and oblate ellipsoids are used as archetypes for
column and plate like ice crystals in clouds ([61]). The in-cloud turbulence
is often not able to destroy the strong alignment of these particles ([9]),
which is in agreement with the orientation model of [34]. As a result, their
orientation statistics and sedimentation velocities can have important
consequences for remote sensing applications like polarization LIDAR, which is
a key component of climate-research programs to characterize the properties of
mixed-phase cloud systems, such as cirrus clouds ([50],[68],[26]). On a side
note, the strong alignment of ice-crystals in the atmosphere can also be
observed by eye since it causes optical phenomena by scattering light, the
origin of the Perry Arc ([67]).

The present study analyses the orientation dynamics of anisotropic particles
settling in a turbulent flow using a combination of analytical, numerical and
experimental approaches. Recently Kramel [36], Menon et al. [45], Menon [46]
have proposed a “rapid-settling theory” for studying orientation dynamics, a
regime wherein the decorrelation time for a Kolmogorov eddy is much larger
than the time a particle takes to settle through a Kolmogorov eddy. Gustavsson
et al. [21] and Anand et al. [1] have also explored the orientation dynamics
of spheroids in the rapid-settling regime using numerical simulations, where
they confirm the transition from random to horizontally aligned orientation
distributions with increasing settling speeds. Gustavsson et al. [21] modelled
the turbulent background flow as a sum of Fourier modes, the kinetic
simulation model. They obtained a normal probability distribution function
(PDF) for the component of the orientation vector along gravity with a
variance that scaled as $S_{v}^{-4}$, in agreement with our earlier studies
[36, 45, 46]. Anand et al. [1] analyzed the sedimentation of spheroids in an
ambient homogeneous isotropic turbulent field, where the background flow was
obtained using direct numerical simulations (DNS). Their calculation of
variance agrees with the previous investigations. Using analyses of the higher
moments, they showed that the orientation PDFs are non-Gaussian due to the
non-Gaussian nature of the turbulent velocity gradient stemming from the
dissipation range intermittency.

In this paper, we propose a new way of investigating non-spherical, inertial
particles which enables us to observe the transition from strongly aligned
particles ($S_{F}\gg 1$) to almost randomly oriented particles ($S_{F}\ll 1$).
Instead of ellipsoidal shaped particles, whose full solid-body rotation is
difficult to measure, we introduce ramified particles. A ramified particle
consists of any number of individual but connected fibers and can be used to
model a wide variety of shapes by adjusting the number and length of the
fibers. A triad, three coplanar fibers connected with equal angles at their
ends, is a crude approximation of a disk-like particle, whereas a jack, three
orthogonal fibers connected at their center, is an approximation for a
spherical-particle. Adjusting the length of each fiber gives us control over
the effective aspect ratio of the ramified particle and in the limit of many
fibers, the ramified particle will approach its ellipsoidal counterpart. The
advantage of using ramified particles over ellipsoidal shaped particles is
that we can measure the orientation very precisely. Moreover, ignoring the
interactions between individual fibers of a ramified particle and using the
well developed slender body theory for single fiber motion, is a good starting
point for more accurate models.

Figure 1: Sedimenting Triad. The particle is leaving a trail of dye behind it.

Section 2 of this paper is a presentation of the theory and stochastic
simulations describing the orientation of thin settling particles in isotropic
turbulence. Starting with a treatment of fibers that are smaller than the
Kolmogorov length scale and settle with small but non-zero Reynolds number, we
progress to consider disk-like particles and triads formed by connecting three
fibers and finally describe modifications of the theory to account for larger
particle sizes and larger particle Reynolds numbers. Section 3 describes the
complementary experimental investigation. First, we discuss the experimental
methods including the turbulent channel apparatus and characterization of the
turbulence as well as the synthesis and characterization of the particles. We
then present results for the orientational behavior of particles settling in
quiescent fluids and in turbulent flows and compares the latter results with
the theoretical predictions. Section 4 is a conclusion and summary of the
study.

##  2 Theory

In this section we present theoretical models for the orientation of
sedimenting fibers and triads in isotropic turbulence. In section 2.1, we
first consider fibers smaller than the Kolmogorov scale that settle with small
but finite $\textit{Re}_{\ell}$. We then outline a method to extend this
theory to larger particle lengths by defining an empirical settling factor
using input from the settling of large particles in a quiescent fluid. The
results in section 2.1 are based on stochastic simulations of the fiber
dynamics in a Lagrangian reference frame. However, fully analytical results
are obtained for fibers settling rapidly through the Kolmogorov scale eddies
corresponding to the limit $S_{F}\gg 1$ in section 2.2. In section 2.3 we
study the rotational dynamics of rapidly settling disks motivated by the
similarity in symmetry of triads and disks. Building upon our understanding of
fibers and disks, a model for the dynamics of small triads is presented in
section 2.4. In section 2.5 we outline approximate modification for the
equations for the translational and rotational motion of fibers and triads for
larger particle Reynolds numbers.

###  2.1 Sedimentation of Small Fibers

In this subsection we present equations governing the orientation and settling
velocity of small, slender fibers. It is assumed that the fiber Reynolds
number $\textit{Re}_{\ell}=W_{min}l/\nu$ is small so that we include inertial
effects only when they break the degeneracy of Stokes flow behavior. Here,
$W_{min}$ is the settling velocity of the particle in a broadside orientation
and $l=L/2$ is the half-length of the particle. The fiber length is much less
than the Kolmogorov length scale $L\ll\eta$, so that the turbulent velocity
field can be approximated as a local linear flow field. In the absence of
particle inertia, fibers experience no net force or torque. In a quiescent
fluid, their orientations, described by a unit vector $\boldsymbol{p}$, are
independent of time and determined by the initial conditions. Fluid inertia
will break this degeneracy and so we include the first effect of fluid inertia
on the hydrodynamic torque experienced by a settling fiber even though
$\textit{Re}_{\ell}$ is small.

Since a settling particle of length $\textit{L}=2\ell$ disturbs a fluid volume
of $\textit{O}(L^{3})$, the particle mass is small compared with the fluid
mass disturbed, if

$\frac{\rho_{p}}{\rho_{f}}\ll\left({\frac{L}{D}}\right)^{2}={\kappa}^{2}$ (1)

From this relation we can see that for high aspect ratio particles, whenever
the particle and fluid densities are comparable, the particle inertia is
negligible compared to that of the fluid as is the case in the experimental
study. We will use this observation to justify the neglect of particle inertia
so that particles experience no net force or torque.

Any external force is balanced by drag and lift forces. Batchelor [3] derived
analytical expressions for the drag and lift force valid for
$\textit{Re}_{\ell}\ll 1$, $\textit{Re}_{D}\ll 1$ and $\kappa\gg 1$. The
balance of forces expressed to leading order in aspect ratio at low
$\textit{Re}_{\ell}$ is given by

$-\frac{4\pi\mu
L}{\ln(2\kappa)}\left(\mathbb{1}-\frac{1}{2}\boldsymbol{p}\boldsymbol{p}\right)\cdot\boldsymbol{W}+m\boldsymbol{g}=0,$
(2)

where $\mathbb{1}$ is the identity matrix, $m=(\rho_{p}-\rho_{f})\pi LD^{2}/4$
is the mass difference between a cylindrical fiber and the displaced fluid,
$\mu$ is the dynamic fluid viscosity and $\boldsymbol{g}$ is the gravitational
acceleration. A fiber will therefore translate with a quasi-steady state
velocity $\boldsymbol{W}$ relative to the local fluid velocity. Equation 2
yields a well-known result for the transverse and longitudinal settling
velocities of a fiber

$W^{f}_{\textit{max}}=2W^{f}_{\textit{min}}$ (3)

where $W^{f}_{\textit{max}}=|\boldsymbol{W}|_{\theta{=}0}$ and
$W^{f}_{\textit{min}}=|\boldsymbol{W}|_{\theta{=}\pi/2}$, respectively. Here,
$\theta$ is the angle between $\boldsymbol{p}$ and $\boldsymbol{g}$.

While one can neglect inertial effects on the settling velocity of small
fibers, fluid inertia will break the degeneracy of particle orientation when a
particle settles in a quiescent fluid. With the inclusion of fluid inertia,
fibers experience inertial torques $\boldsymbol{G}_{sed}$ that rotate the
particle to an equilibrium orientation where $\boldsymbol{p}$ is perpendicular
to $\boldsymbol{W}$. Khayat and Cox [31] derived expressions for the torque
experienced by a translating fiber $\boldsymbol{G}_{sed}$ (see their Eq.
6.12), which becomes in the low Reynolds number limit
$\left(\textit{Re}_{\ell}\ll 1\right)$ and to leading order in small aspect
ratio

$\boldsymbol{G}_{sed}=\frac{5\pi\rho_{f}L^{3}}{24(\ln
2\kappa)^{2}}(\boldsymbol{W}\cdot\boldsymbol{p})(\boldsymbol{W}\times\boldsymbol{p})$
(4)

The particle also experiences a rotational resistance $\boldsymbol{G}_{rel}$
to its relative rotation [3]:

$\displaystyle\boldsymbol{G}_{\textit{rel}}=-\frac{\pi\mu
L^{3}}{3\ln(2\kappa)}(\mathbb{1}-\boldsymbol{p}\boldsymbol{p}).\boldsymbol{\Omega}_{\textit{rel}}$
(5)

Here, $\boldsymbol{\Omega}_{\textit{rel}}$ is the rotation of the particle
relative to the local fluid rotation. In addition, fibers experience torques
due to the fluid strain rate
$\boldsymbol{S}=\frac{1}{2}(\boldsymbol{\Gamma}+\boldsymbol{\Gamma}^{T})$:

$\boldsymbol{G}_{strain}=\frac{\pi\mu
L^{3}}{3\ln(2\kappa)}(\boldsymbol{p}\times(\boldsymbol{S}\cdot\boldsymbol{p}))$
(6)

Here, $\Gamma_{ij}=\partial u_{i}/\partial x_{j}$ is the turbulent velocity
gradient. For a symmetric fiber sedimenting in turbulence, in the absence of
particle inertia, a torque balance reads:

$\underbrace{\frac{5\pi\rho_{f}L^{3}}{24(\ln
2\kappa)^{2}}(\boldsymbol{W}\cdot\boldsymbol{p})(\boldsymbol{W}\times\boldsymbol{p})}_{\textit{inertial
sedimentation}}-\underbrace{\frac{\pi\mu
L^{3}}{3\ln(2\kappa)}(\mathbb{1}-\boldsymbol{p}\boldsymbol{p}).\boldsymbol{\Omega}_{\textit{rel}}}_{\textit{relative
rotation}}+\underbrace{\frac{\pi\mu
L^{3}}{3\ln(2\kappa)}(\boldsymbol{p}\times(\boldsymbol{S}\cdot\boldsymbol{p}))}_{\textit{turbulent
strain}}=0.$ (7)

The torque balance yields the following equation for the time rate of change
$\boldsymbol{\dot{p}}$ of the fiber orientation

$\boldsymbol{\dot{p}}=\boldsymbol{\Gamma}.\boldsymbol{p}-\boldsymbol{p}\left(\boldsymbol{p}.\boldsymbol{S}.\boldsymbol{p}\right)+\frac{5}{8\nu\ln{(2\kappa)}}(\boldsymbol{W}{\cdot}\boldsymbol{p})~{}\boldsymbol{W}{\cdot}(\boldsymbol{p}\boldsymbol{p}-\mathbb{1})$
(8)

where the first two terms correspond to Jeffery rotation in the local linear
flow field [29] and the last term is the rotation due to the inertial torque
caused by the particles sedimentation. Without a background flow, the inertial
torque acts to orient a sedimenting spheroidal particle broadside-on to
gravity. However, in the presence of additional torque, such as a
gravitational torque for an axisymmetric particle with mass asymmetry, the
inertial torque can compete to create an oblique settling orientation [55]. In
the current study, the torque due to turbulent shear competes with the
inertial torque.

A time scale $\tau_{sed}$ of rotation due to sedimentation of a fiber at small
$\textit{Re}_{\ell}$ may be defined as the inverse rotation rate of the
particle in quiescent fluid at $\theta=45^{\circ}$, where $\theta$ is the
angle between $\boldsymbol{p}$ and $\boldsymbol{g}$. When we generalize this
definition to disk and triad particles, $\boldsymbol{p}$ will be the normal
vector to the plane of the particle. Upon solving Eq. 2 and the component of
the tumbling rate due to the inertial torque (from Eq. 8), at
$\theta=45^{\circ}$ we get,

$\tau_{\textit{sed}}=\frac{8\nu\ln\left(2\kappa\right)}{5W^{2}_{min}}$ (9)

where we have used Eq. 2 to find $W_{min}$, the minimum settling velocity of a
fiber which is achieved at $\theta=90^{\circ}$. The typical timescale of
response of small fibers due to turbulence $\tau_{\textit{turb}}$ is the
Kolmogorov timescale,

$\tau_{\eta}=\frac{\eta}{u_{\eta}}=\sqrt{\frac{\nu}{\epsilon}}$ (10)

The orientation distribution of fibers may now be understood by comparing the
two time scales $\tau_{\eta}$ and $\tau_{\textit{sed}}$. We define the
settling factor as the ratio of these two time scales to be,

$S_{F}=\frac{\tau_{\eta}}{\tau_{\textit{sed}}}$ (11)

When $S_{F}\gg 1$, the rotation due to turbulence is weak compared to that due
to the inertial torque and leads to a small deviation from the horizontal
orientation, making $\boldsymbol{W}$ align parallel to gravity and
perpendicular to $\boldsymbol{p}$. On the other hand at $S_{F}\ll 1$,
turbulence is relatively stronger leading to an isotropic orientation
distribution. For small fibers, ($L\ll\eta$), and $\textit{Re}_{\ell}\ll 1$
the expression for $S_{F}^{f}$ reduces to,

$S^{f}_{F}=\frac{5W^{2}_{min}\tau_{\eta}}{8\nu\ln{(2\kappa)}}=\frac{5}{8\ln{(2\kappa)}}\left(\frac{W_{min}}{u_{\eta}}\right)^{2}$
(12)

The superscript $f$ in Eq. 12 indicates the settling factor for small fibers.
We will continue to use the general definition Eq. 11 to define settling
factors for small triads and disks in the subsequent developments. In the
cloud microphysics literature, the parameter Sv denotes the ratio of the
Kolmogorov eddy turnover time to the time a particle takes to settle across an
eddy, which can also be written as $S_{v}=W/u_{\eta}$ [12]. Thus, we can see
that the settling factor at small Reynolds number is proportional to the
square of this non-dimensional settling velocity, i.e., $S_{F}^{f}\propto
S_{v}^{2}$.

To determine the fiber orientation, we must solve Eq. 2 and Eq. 8 in a
reference frame translating with the particle. For slowly settling fibers
$S_{F}\ll 1$, the fiber follows a Lagrangian path. For rapidly settling fibers
$S_{F}\gg 1$, it will be shown in the next subsection that the fiber
orientation arises from a quasi-steady balance of turbulent shear and inertial
rotation so that the particle path does not change the orientation. One might
then expect to obtain a reasonable estimate of the fiber orientation by
simulating fibers experiencing the turbulence on a Lagrangian path for all
$S_{F}$. For this purpose we employ a stochastic model. Meneveau [44] reviews
stochastic models to describe the fluid velocity gradient along a Lagrangian
path. We have employed the model developed by Girimaji and Pope [19] that
captures the log-normal distribution of the pseudo-dissipation, the time scale
for relaxation of the strain rate tensor on a Lagrangian path, and the
tendency of the nonlinear inertial terms to align the vorticity with the
strain axes. Shin and Koch [57] showed that this model predicts the rotational
velocity variance of neutrally buoyant particles computed in direct-numerical
simulations (DNS) with much greater accuracy than a simple Gaussian velocity
gradient model (Brunk et al. [7]). Girimaji and Pope obtained favorable
comparisons of many of the tensor invariants of turbulence with the DNS of
Yeung and Pope [71]. The inputs to the model which include the correlation
times for the pseudo-dissipation and the components of the strain rate and the
variance of the logarithm of the pseudo-dissipation can be obtained from Yeung
and Pope for several values of $Re_{\lambda}$, and we use $Re_{\lambda}=38$
and $93$ for our simulations. Recently, in another problem of particles in a
turbulent flow, we have used the Lagrangian velocity gradient model of
Girimaji and Pope [19] to explore the role of non-Gaussian statistics in
collisions of hydrodynamically interacting particles settling in a turbulent
flow [13].

In Figure 2, the orientational variance for a fiber settling in a turbulent
flow with $Re_{\lambda}=38$ is shown as a function of the settling factor. It
is seen that $\langle\cos^{2}(\theta)\rangle$, where $\theta$ is the angle
between $\boldsymbol{p}$ and $\boldsymbol{g}$, shows a smooth transition from
an isotropic orientation distribution
$\langle\cos^{2}(\theta)\rangle=\frac{1}{3}$ to nearly aligned distribution in
which orientational dispersion about the horizontal plane ($\cos\theta=0$)
decays as a power law in the settling factor. As will be derived analytically
in section 2.2, the decay in the orientational variance about the equilibrium
orientation $\theta=90^{\circ}$ follows,

$\langle\cos^{2}\theta\rangle=0.022~{}{S^{f}_{F}}^{-2}$ (13)

The scaling of the orientational variance with $S_{F}^{f}$ is in agreement
with the findings of Gustavsson et al. [21]
($\langle\cos^{2}\theta\rangle\propto$ Sv-4) and Anand et al. [1]. Anand et
al. [1] evaluates the variance in terms of a Froude number, Frη, that is
identical in definition to that of Sv.

The theory can be extended in an approximate way to larger particles by
replacing the time scale of rotation due to sedimentation
$\tau_{\textit{sed}}$ with $T_{\textit{sed}}$ extracted from experimental
observations of a particle settling in a quiescent fluid

$T_{\textit{sed}}=\frac{1}{|\boldsymbol{\dot{p}}|}_{\theta=45^{\circ}}$ (14)

Moreover, the rotations of particles larger than the Kolmogorov length scale
are dominated by turbulent eddies close to their size. Therefore, the
appropriate turbulent time scale for rotations is the eddy turn over time at
the scale of the particle (Parsa and Voth [53])

$\tau_{L}=\frac{L}{u_{L}}=\sqrt{\frac{4}{15}}\frac{L}{u_{L}^{T}},$ (15)

where
$u_{L}^{T}=\sqrt{\langle(\Delta\mathbf{u}\cdot(\mathbb{1}-\boldsymbol{\hat{r}}\boldsymbol{\hat{r}}))^{2}\rangle}$
is the root-mean-square of the transverse components of the fluid velocity
difference
$\langle(\Delta\mathbf{u}=\mathbf{u}(\mathbf{x}+L\boldsymbol{\hat{r}})-\mathbf{u}(\mathbf{x})$)
at the scale of the particle and the factor of $\sqrt{4/15}$ is chosen such
that $\tau_{L}\to\tau_{\eta}$ in the limit $L\ll\eta$. The ratio of these two
time scales allows us to define a settling factor

$S_{F}=\frac{\tau_{L}}{T_{sed}}$ (16)

that will be used characterize relative rates of sedimentation-induced and
turbulent-induced rotation rate in the experiments.

Figure 2: Orientation variance of small fibers as a function of the settling
factor, $S_{F}$. The squares correspond to simulations and the lines are
asymptotes derived in the low and high $S_{F}$ limits. The solid line at high
$S_{F}$ is the rapid settling limit in a particle reference frame for which
the particle orientation and velocity gradient are uncorrelated. The dashed
line is the asymptote for a Lagrangian frame which captures the correlation of
transverse particle orientation with velocity gradient observed in the
simulation.

###  2.2 Rapid Settling of Small Fibers

In this subsection, we will derive an analytical prediction for the variance
of the orientation in the rapid settling limit,
$S_{F}=\tau_{\eta}/\tau_{\textit{sed}}\gg 1$. This indicates that the
relaxation of fiber orientation toward its equilibrium horizontal alignment
occurs much more rapidly than the Kolmogorov time scale. From Eq. 12, it can
be seen that this limit corresponds to one in which the time
$\tau_{\textit{samp}}=\eta/W$ for a fiber to sample a Kolmogorov scale eddy is
also much smaller than $\tau_{\eta}$. In particular,
$\tau_{\eta}/\tau_{\textit{samp}}=W/u_{\eta}={S^{f}_{F}}^{1/2}[\ln(2\kappa)]^{1/2}\gg
1$. From these results it can be seen that the relaxation of fiber orientation
is much faster than the sampling time

$\frac{\tau_{\textit{sed}}}{\tau_{\textit{samp}}}=\frac{[\ln(2\kappa)]^{1/2}}{{S^{f}_{F}}^{1/2}}\ll
1$ (17)

Thus, despite the rapid translational motion of the particle through the
eddies, the fiber responds to changes in the local shear rate and achieves a
new orientation sufficiently rapidly so that one may obtain the orientation
from a quasi-steady balance of the rotation due to sedimentation and turbulent
shear.

We will determine the rotation rate of fibers whose orientations exhibit small
deviations from the horizontal plane, so that $\langle{p_{3}}^{2}\rangle\ll 1$
where the $3$ axis is parallel to gravity. We begin with an alternate mobility
form of Eq. 2 written using Einstein notation as

$W_{i}=\frac{\ln{2\kappa}}{8\pi\mu
l}\left(\delta_{ij}+p_{i}p_{j}\right)F\delta_{j3}$ (18)

where $F=mg$. For $\langle{p_{3}}^{2}\rangle\ll 1$, it can be seen that

$W_{i}p_{i}=2W_{3}p_{3}=\frac{\ln{2\kappa}}{8\pi\mu l}Fp_{3}$ (19)

Substituting this result into Eq. 8 yields,

$\dot{p}_{i}=\frac{10Re_{\ell}p_{3}}{8\ell\ln{2\kappa}}2W_{3}p_{3}p_{i}-\frac{10Re_{\ell}p_{3}}{8\ell\ln{2\kappa}}W_{i}+\Gamma_{ij}p_{j}-p_{i}S_{jl}p_{j}p_{l}$
(20)

Since the bodies will remain horizontal on average in this limit, as expected
from theory and shown in simulation, we have $p_{3}\ll p_{1,2}$. Thus,

$\dot{p}^{t}_{i}=\Gamma_{ij}p^{t}_{j}-p^{t}_{i}S_{jl}p^{t}_{j}p^{t}_{l}$ (21)

where $p^{t}_{i}$ denotes the transverse component ($i=1,2$) of the
orientation vector, indicating that the fiber orientation samples the plane
normal to gravity by turbulent shearing motions. This will lead to an
isotropic distribution of orientation within the 1-2 plane. Since
$\tau_{\textit{sed}}\ll\tau_{\eta}$, the motion within the 1-2 plane will be
slow compared with the equilibration of the 3 component of the fiber
orientation with the current turbulent shear flow. In the settling direction
(3) there will exist a quasi-static equilibrium because
$\tau_{sed}\ll\tau_{\textit{samp}}$. Thus,

$\dot{p}_{3}=\frac{10Re_{\ell}p_{3}}{8\ell\ln{2\kappa}}2W_{3}p_{3}p_{3}-\frac{10Re_{\ell}p_{3}}{8\ell\ln{2\kappa}}W_{3}+\delta_{i3}\Gamma_{ij}p^{t}_{j}-p_{3}S_{jl}p^{t}_{j}p^{t}_{l}=0$
(22)

Since $p_{3}\ll 1$ we can further simplify the above equation by balancing the
second and third terms to obtain,

$p_{3}\sim\frac{8\ell\ln{2\kappa}}{10Re_{\ell}W_{3}}\delta_{i3}\Gamma_{ij}p^{t}_{j}$
(23)

Thus, the variance characterizing the ”wiggle” out of the horizontal plane is

$\langle
p_{3}^{2}\rangle=\left(\frac{8\ell\ln{2\kappa}}{10Re_{\ell}W_{3}}\right)^{2}\delta_{i3}\delta_{m3}\langle
p^{t}_{j}p^{t}_{n}\rangle\langle\Gamma_{ij}\Gamma_{mn}\rangle=\left(\frac{8\ell\ln{2\kappa}}{10Re_{\ell}W_{3}}\right)^{2}\delta_{i3}\delta_{m3}\langle
p^{t}_{j}p^{t}_{n}\rangle\left[\langle S_{ij}S_{mn}\rangle+\langle
R_{ij}R_{mn}\rangle\right]$ (24)

where the $\langle\rangle$ denotes ensemble averages and
$R_{ij}=1/2(\Gamma_{ij}-\Gamma_{ji})$ is the antisymmetric part of the
velocity gradient. Cross terms such as $\langle S_{ij}R_{mn}\rangle$ are zero
due to isotropy of the turbulent field. In the above expression, we have
assumed $\langle p^{t}_{j}p^{t}_{n}\Gamma_{3j}\Gamma_{3n}\rangle$ to be the
corresponding product of mean of velocity gradient
$\langle\Gamma_{3j}\Gamma_{3n}\rangle$ and orientation $\langle
p^{t}_{j}p^{t}_{n}\rangle$ dyads. This is a valid assumption in the rapid
settling limit because $\tau_{sed}\ll\tau_{\eta}$ and, as a result,
$\boldsymbol{p}^{t}$ changes on a much larger time scale than the velocity
gradient. $\langle S_{ij}S_{mn}\rangle$ and $\langle R_{ij}R_{mn}\rangle$ are
fourth order isotropic tensors whose form can be deduced using the properties
of symmetry and continuity, as shown in Brunk et al. [7] to obtain

$\langle
S_{ij}S_{mn}\rangle=\frac{S^{2}}{10}\left[\delta_{im}\delta_{jn}+\delta_{in}\delta_{jm}-\frac{2}{3}\delta_{ij}\delta_{mn}\right]$
(25a) $\langle
R_{ij}R_{mn}\rangle=\frac{R^{2}}{6}\left[\delta_{im}\delta_{jn}-\delta_{in}\delta_{jm}\right]$
(25b)

where $S^{2}=\langle S_{ij}S_{ij}\rangle$ and $R^{2}=\langle
R_{ij}R_{ij}\rangle$. Substituting Eq. 25b in Eq. 24 we have

$\displaystyle\langle p_{3}^{2}\rangle$
$\displaystyle=\left(\frac{8\ell\ln{2\kappa}}{10Re_{\ell}W_{3}}\right)^{2}\delta_{i3}\delta_{m3}\left[\frac{S^{2}}{10}\left(\delta_{im}+\frac{\langle
p_{i}p_{m}\rangle}{3}\right)+\frac{R^{2}}{6}\left(\delta_{im}-\langle
p_{i}p_{m}\rangle\right)\right]$
$\displaystyle=\left(\frac{8\ell\ln{2\kappa}}{10Re_{\ell}W_{3}}\right)^{2}\left[\frac{S^{2}}{10}+\frac{R^{2}}{6}\right]$
$\displaystyle\langle
p_{3}^{2}\rangle=\left(\frac{8\ell\ln{2\kappa}}{5Re_{\ell}W_{3}}\right)^{2}\frac{\Gamma_{\eta}^{2}}{30}=\frac{1}{30}{S^{f}_{F}}^{-2}$
(26)

where we have used the relations $S^{2}=R^{2}=\Gamma_{\eta}^{2}/2$ for
homogeneous isotropic turbulence and $S^{f}_{F}$ from Eq. 12. Thus, we have
for the rapid settling limit the following relation characterizing the
departure of orientation from the horizontal plane due to turbulence,

$\langle\cos^{2}\theta\rangle=\frac{1}{30}{S^{f}_{F}}^{-2}\approx
0.033~{}{S^{f}_{F}}^{-2}$ (27)

In our simulations we use a Langrangian model of velocity gradient instead of
a particle frame model of turbulence. In the large $S^{f}_{F}$ limit of the
simulations, while the strong inertial torque tries to maintain a horizontal
orientation, the fiber orientation continues to change on the Kolmogorov time
scale in the 1-2 plane. This leads to a correlation of the transverse fiber
orientation with the velocity gradient. As a result, our simulations are
different from Eq. 27 by a factor corresponding to $\langle
p^{t}_{i}p^{t}_{j}\Gamma_{3i}\Gamma_{3j}\rangle/\left(2\Gamma_{\eta}^{2}/15\right)$.
In our simulations, this factor is observed to be around 0.66, making the
asymptote

$\langle\cos^{2}\theta\rangle=\frac{0.66}{30}{S^{f}_{F}}^{-2}\approx
0.022{S^{f}_{F}}^{-2}$ (28)

###  2.3 Rapid settling of Disks

In this subsection, we will derive an analytical prediction for the variance
of the orientation of rapidly settling disks before studying the rotational
dynamics of triads in 2.4. Disks and triads are geometrically similar - triads
have 3-fold rotational symmetry while disks have circular symmetry. This hints
at possible similarities of the resistance tensors of the two objects. (See
Brenner [6] for a general discussion on the equality of resistance tensors for
objects based on symmetry considerations.)

An oblate spheroid of semi-major axis length $l$ spinning with relative
angular velocity $\bm{\Omega}_{rel}$ with respect to the local fluid rotation
in a linear flow experiences a net torque

$\displaystyle\mathbf{G}_{rel+strain}=-8\pi\mu
l^{3}\left[X^{C}\mathbf{pp}+Y^{C}(\mathbb{1}-\mathbf{pp})\right].\bm{\Omega}_{rel}+8\pi\mu
l^{3}Y^{H}(\mathbf{p}\times(\mathbf{S.p})).$ (29)

Here $X^{C},Y^{C}$ and $Y^{H}$ are the scalar resistance functions associated
with the dynamics of a spheroidal particle in Stokes flow [32]. For thin disks
they take the values

$\displaystyle
X^{C}\sim\frac{4}{3\pi},\,\,Y^{C}\sim\frac{4}{3\pi},\,\,Y^{H}\sim-\frac{4}{3\pi}$
(30)

With the inclusion of fluid inertia, disks, like fibers, experience inertial
torques that rotate toward an equilibrium orientation with the large
dimensions of the particle perpendicular to the velocity. In the case of
disks, this corresponds to $\mathbf{p}$ aligned with $\mathbf{W}$. Dabade et
al. [11] derived expressions for the torque experienced by a sedimenting
spheroidal particle, assuming fluid inertia to be weak. For oblate spheroids,
they derive the torque on a thin disk of semi-major axis length $l$

$\displaystyle\mathbf{G}_{sed}=-\left\\{\frac{38}{9}-\frac{17216}{945\pi^{2}}\right\\}\frac{l^{3}\rho_{f}}{8}(\mathbf{W.p})(\mathbf{W}\times\mathbf{p})$
(31)

For a thin disk sedimenting in turbulence, in the absence of particle inertia,
a torque balance reads:

$\displaystyle-\underbrace{\left\\{\frac{38}{9}-\frac{17216}{945\pi^{2}}\right\\}\frac{l^{3}\rho_{f}}{8}(\mathbf{W.p})(\mathbf{W}\times\mathbf{p})}_{\textit{inertial
sedimentation}}-\underbrace{\frac{32\mu
l^{3}}{3}\bm{\Omega}_{rel}}_{\textit{relative
rotation}}-\underbrace{\frac{32\mu
l^{3}}{3}(\mathbf{p}\times(\mathbf{S.p}))}_{\textit{turbulent strain}}=0$ (32)

The zero torque balance yields the following equation for the time rate of
change $\dot{\mathbf{p}}$ of the disk orientation

$\displaystyle\dot{\mathbf{p}}=\mathbf{R}.\mathbf{p}-\mathbf{S}.\mathbf{p}+\mathbf{p}(\mathbf{p}.\mathbf{S}.\mathbf{p})-\frac{c}{\nu}(\mathbf{W}.\mathbf{p})\mathbf{W}.(\mathbf{pp}-\mathbb{1})$
(33)

where $c=(19/32-269/105\pi^{2})/12\approx 0.028$. The first three terms on the
right-hand side of Equation 33 are the Jeffery rotation rate [29] for a
particle with $\kappa\ll 1$ and the third term is the rotation due to the
settling-induced inertial torque.

We will now determine the rotation rate of rapidly settling disks that remain
nearly aligned with the 3 axis with small fluctuations in the horizontal
plane, so that $<p_{3}^{2}>\sim 1$. We begin with the mobility expression
giving the sedimentation velocity

$\displaystyle W_{i}=\frac{3}{32\mu
l}(\delta_{ij}-\frac{1}{3}p_{i}p_{j})F\delta_{j3}$ (34)

where $F=mg$. Equation 34 yields the following well-known result for the
transverse and longitudinal settling velocities of a thin disk

$\displaystyle W^{d}_{max}=1.5W^{d}_{min}=3F/(32\mu l)$ (35)

where $W^{d}_{max}=|\mathbf{W}|_{\theta=\pi/2}$ and
$W^{d}_{min}=|\mathbf{W}|_{\theta=0}$, respectively. For $<p_{3}^{2}>\sim 1$

$\displaystyle W_{i}p_{i}=W_{3}p_{3}=\frac{Fp_{3}}{16\mu l}$ (36)

Substituting this result into Eq. 33 yields

$\displaystyle\dot{p}_{i}=R_{ij}p_{j}-S_{ij}p_{j}+p_{i}p_{k}S_{kl}p_{l}+\frac{3cW_{3}^{2}}{2\nu}(\delta_{ij}-p_{i}p_{j})\left(\delta_{j3}-\frac{1}{3}p_{j}p_{3}\right)p_{3}$
(37)

Since the disk remains nearly horizontal, we can write
$p_{i}=\delta_{i3}+p_{i}^{t}$ where $p_{i}^{t}\ll 1$. The rotation rate of the
transverse component of the orientation vector, $p^{t}_{i}$ ($i=1,2$), then
assumes the following simplified form

$\displaystyle\dot{p}^{t}_{i}\approx
R_{i3}-S_{i3}-\frac{3cW_{3}^{2}}{2\nu}p^{t}_{i}=0$ (38)
$\displaystyle\Rightarrow$ $\displaystyle
p^{t}_{i}\approx\frac{2\nu}{3cW_{3}^{2}}(R_{i3}-S_{i3}).$ (39)

Thus the variance characterizing the “wiggle” out of the vertical axis is

$\displaystyle\left\langle p^{t\,2}_{i}\right\rangle$ $\displaystyle=$
$\displaystyle\left(\frac{2\nu}{3cW_{3}^{2}}\right)^{2}\left(\left\langle
R_{i3}R_{i3}\right\rangle+\left\langle
S_{i3}S_{i3}\right\rangle\right)=\left(\frac{2\nu}{3cW_{3}^{2}}\right)^{2}\left(\frac{\Gamma^{2}_{\eta}}{6}+\frac{\Gamma^{2}_{\eta}}{10}\right)$
(40) $\displaystyle=$
$\displaystyle\left(\frac{2\nu}{3cW_{3}^{2}}\right)^{2}\frac{\Gamma^{2}_{\eta}}{3}=\left(\frac{2\nu}{3cW^{d\,2}_{min}}\right)^{2}\frac{4\Gamma^{2}_{\eta}}{15}$

Similar to fibers we can define a settling factor for disks as

$\displaystyle
S_{F}^{d}=\frac{\tau_{\eta}}{\tau_{sed}}=\frac{3cW^{d\,2}_{min}}{4\Gamma_{\eta}}$
(41)

where $\tau_{sed}=4/(3cW^{d\,2}_{min})$ is the time scale of rotation due to
sedimentation for a disk whose axis is at a 45∘ angle to gravity. The small
orientation variance in the high $S^{d}_{F}$ limit is

$\displaystyle 1-\left\langle
p_{3}^{2}\right\rangle=\left\langle\sin^{2}\theta\right\rangle=\left\langle
p^{t\,2}_{i}\right\rangle=\frac{1}{15}S_{F}^{d\,-2}$ (42)

To compare orientational behavior of disks with fibers, we define the average
deviation of a disk away from horizontal, similar to the
$\langle\cos^{2}\theta\rangle$ for fibers,

$0.50\left(1-\langle p_{3}^{2}\rangle\right)=\frac{1}{30}{S^{d}_{F}}^{-2}$
(43)

###  2.4 Triads Settling in Turbulence

In this subsection, the theory and simulation are extended to triads - three
armed ramified particles where all three fiber arms lie in the same plane at
equal separation. We model ramified particles as hydrodynamically independent
fibers connected together to translate and rotate as a single rigid body.
Thus, the force and torque balances on the triad are:

$\sum_{n=1}^{3}[\boldsymbol{F}^{n}_{drag}+\boldsymbol{F}^{n}_{gravity}]=\sum_{n=1}^{3}\left[-\frac{4\pi\mu
L}{\ln(2\kappa)}\left(\mathbb{1}-\frac{1}{2}\boldsymbol{p}^{\prime}_{n}\boldsymbol{p}^{\prime}_{n}\right)\cdot\boldsymbol{W}_{n}+m_{n}\boldsymbol{g}\right]=0$
(44a) $\displaystyle\sum_{n=1}^{3}$
$\displaystyle\left[\underbrace{\frac{5\pi\rho_{f}L^{3}}{24(\ln
2\kappa)^{2}}(\boldsymbol{W}_{n}{\cdot}\boldsymbol{p}^{\prime}_{n})(\boldsymbol{W}_{n}{\times}\boldsymbol{p}^{\prime}_{n})}_{inertial\,sedimentation}-\underbrace{\frac{\pi\mu
L^{3}}{3\ln(2\kappa)}(\mathbb{1}{-}\boldsymbol{p}^{\prime}_{n}\boldsymbol{p}^{\prime}_{n}){\cdot}\boldsymbol{\Omega}_{\textit{rel}}}_{relative\,rotation}\right.$
$\displaystyle+\left.\underbrace{\frac{\pi\mu
L^{3}}{3\ln(2\kappa)}(\boldsymbol{p}^{\prime}_{n}{\times}(\boldsymbol{S}\cdot\boldsymbol{p}^{\prime}_{n}))}_{turbulent\,strain}+\underbrace{\ell\boldsymbol{p}^{\prime}_{n}\times\boldsymbol{F}^{n}_{drag}}_{drag\,on\,arms}\right]=0$

where $L=2\ell$ is the arm length,$\boldsymbol{W}^{c}$ is the relative
velocity of the triad center of mass with the fluid and
$\boldsymbol{W}_{n}=\boldsymbol{W}^{c}+\ell\boldsymbol{\Omega}^{c}\times\boldsymbol{p}_{n}^{\prime}-\ell\boldsymbol{p}_{n}^{\prime}.\boldsymbol{\Gamma}$
is the relative velocity of the nth arm with the fluid. The orientation of
each arm is defined by $\boldsymbol{p}_{n}^{\prime}$ and the orientation of a
ramified particle is defined by $\boldsymbol{p}$, which is perpendicular to
the plane spanned by the arms. From the symmetry of the particle, the
gravitational torque sums to zero. As in the case of disks, the minimum
velocity occurs when $\boldsymbol{p}$ is parallel to gravity, so that
$W^{t}_{\textit{max}}=|\boldsymbol{W}|_{\theta{=}\pi/2}$ and
$W^{t}_{\textit{min}}=|\boldsymbol{W}|_{\theta{=}0}$.

This model neglects hydrodynamic interactions among the rods. The influence of
hydrodynamic interactions on the drag and the torques due to relative rotation
and straining motions are of higher order in the small parameter
$1/\ln(2\kappa)$. However, hydrodynamic interactions would influence the
inertial torque at the same order of magnitude as the terms retained and the
present model that includes only the torque on each arm acting independently
is likely an underestimate of the triads’ inertial torque. When comparing with
experimental measurements of orientation in turbulent flows, the
experimentally observed inertial rotation rate of a large triad is used to
correct for this discrepancy.

It is important to note that our theory applies to a case where the Reynolds
number is small $Re_{\ell}\ll 1$. In this limit, the rotation of the triad
toward horizontal orientations is slow. The competition between turbulent
shear and inertial rotation leads to intermediate orientation distributions
between isotropic and full alignment when the turbulence is weak
$G=\ell\Gamma_{\eta}/W^{f}_{min}\ll 1$ and the Reynolds number is small
$Re_{\ell}\ll 1$, but the settling parameter
$S^{f}_{F}=\left(\frac{5}{16}\right)\frac{Re_{\ell}}{G}=O(1)$. In this limit
the velocity of the triad center of mass $\boldsymbol{W}^{c}$ is much larger
than the relative velocity of the arms with respect to the triad center of
mass, so that the inertial torque due to the translational motion of the
particle dominates that due to the triad’s rotation. In order to ensure these
conditions in our simulations, especially at higher settling rates, we scale
our force and torque balance equations. We scale Eq. 44a using
${\mu}W^{f}_{min}\ell^{2}$ and Eq. LABEL:eq:torque-balance-ram using
${\mu}{\ell^{3}}\Gamma_{\eta}$, and express Eq. 44 with
$\boldsymbol{W}_{n}=\boldsymbol{W}^{c}+\boldsymbol{w}_{n}$, where
$\boldsymbol{W}^{c}$ is the velocity of triad center of mass and
$\boldsymbol{w}_{n}\propto\ell\Gamma_{\eta}$ is the disturbance velocity
experienced by arms $n=1,2,3$. In the limit of $Re_{\ell}\ll 1$ and
$G=\ell\Gamma_{\eta}/W^{f}_{min}\ll 1$, the triad equations reduce to

$\sum_{n=1}^{3}\left[\left(\mathbb{1}-\frac{1}{2}\boldsymbol{p}^{\prime}_{n}\boldsymbol{p}^{\prime}_{n}\right)\cdot\bar{\boldsymbol{W}}^{c}-\hat{\boldsymbol{e}}_{g}\right]=0$
(45a)
$\displaystyle\sum_{n=1}^{3}\left[S^{f}_{F}\left(\bar{\boldsymbol{W}}^{c}{\cdot}\boldsymbol{p}^{\prime}_{n}\right)\left(\bar{\boldsymbol{W}}^{c}{\times}\boldsymbol{p}^{\prime}_{n}\right)-4(\mathbb{1}{-}\boldsymbol{p}^{\prime}_{n}\boldsymbol{p}^{\prime}_{n}){\cdot}\bar{\boldsymbol{\Omega}^{\textit{c}}}\right.$
$\displaystyle+\left.4(\boldsymbol{p}^{\prime}_{n}{\times}(\bar{\boldsymbol{\Gamma}}\cdot\boldsymbol{p}^{\prime}_{n}))\right]=0$
(45b)

In Eq. 45a, $W^{f}_{min}=\frac{mg\ln 2\kappa}{4\pi\mu L}$ is transverse
velocity of a settling fiber,
$\bar{\boldsymbol{W}}^{c}=\frac{\boldsymbol{W}^{c}}{W^{f}_{min}}$, and
$\hat{\boldsymbol{e}}_{g}$ is the gravitational unit vector. In Eq. 45b, the
settling factor is defined using the definition for fibers in Eq. 12,
$\bar{\boldsymbol{\Omega}}^{c}=\frac{\boldsymbol{\Omega}^{c}}{\Gamma_{\eta}}$,
and $\bar{\boldsymbol{S}}=\frac{\boldsymbol{S}}{\Gamma_{\eta}}$. The symmetry
of the triad leads to
$\sum_{n=1}^{3}\boldsymbol{p}^{\prime}_{n}\boldsymbol{p}^{\prime}_{n}=3/2(\mathbb{1}-\boldsymbol{pp})$,
the signature of a body whose hydrodynamic response is transversely isotropic.
The triad equations (Eq. 45a-45b) then assume the following simplified forms -

$\displaystyle\bar{\boldsymbol{W}}^{c}=\frac{4}{3}\left(\mathbb{1}-\frac{1}{4}\boldsymbol{p}\boldsymbol{p}\right).\hat{\boldsymbol{e}}_{g}$
(46)
$\displaystyle-S^{f}_{F}\left(\bar{\boldsymbol{W}}^{c}{\cdot}\boldsymbol{p}\right)\left(\bar{\boldsymbol{W}}^{c}{\times}\boldsymbol{p}\right)-4(\mathbb{1}{-}\boldsymbol{p}\boldsymbol{p}){\cdot}\bar{\boldsymbol{\Omega}^{\textit{c}}}+4(\boldsymbol{p}{\times}(\bar{\boldsymbol{\Gamma}}\cdot\boldsymbol{p}))=0$
(47)

In the low Reynolds number limit, this ramified particle model (Eq. 46)
predicts a different ratio of maximum and minimum sedimentation velocities for
triads than for disks. In place of Eqn. 35, we have

$\displaystyle W^{t}_{\textit{max}}=\frac{4}{3}W^{t}_{\textit{min}}$ (48)

In the rapid settling limit the triad is approximated to be in a quasi-steady
nearly horizontal orientation, allowing the angular velocity that would rotate
the triad out of 12-plane to be neglected. Eq. (47) then simplifies to

$\displaystyle S_{F}^{f}(-\delta_{i1}p_{2}+\delta_{i2}p_{1})\approx
3\epsilon_{imk}\Gamma_{km}-3\epsilon_{ijk}\Gamma_{km}p_{j}p_{m}.$ (49)

Similar to disks, the variance characterizing the “wiggle” of a triad out of
the vertical axis is

$\displaystyle\langle
p^{t\,2}_{i}\rangle=\frac{9}{S_{F}^{f}}\left[\langle\Gamma_{31}^{2}\rangle+\langle\Gamma_{32}^{2}\rangle\right]$
(50)

where $p^{t}_{i}$ ($i=1,2$) is the transverse component of the orientation
vector. The asymptotic expression for a triad is qualitatively similar to that
for a disk, with a power law dependence.

$1-<p_{3}^{2}>=\langle\sin^{2}\theta\rangle=\frac{12}{5}{S^{f}_{F}}^{-2}.$
(51)

However, we see a difference of coefficient compared to Eq. (42) for the disk
orientational moment. This reflects differences in the inertial rotation rate
of triads and disks as well as the fact that we have used $S^{f}_{F}$ for the
settling factor in the preceding development. To define a settling factor for
triads, we use Eq. (47) to obtain an expression for ${\boldsymbol{\dot{p}}}$
in a quiescent fluid and find the rotational time scale for a triad oriented
at a 45∘ angle to gravity to be $\tau_{sed}=6/(S_{F}^{f}\Gamma_{\eta})$. Thus,
the triad settling factor is

$\displaystyle S_{F}^{t}=\frac{\tau_{\eta}}{\tau_{sed}}=\frac{S^{f}_{F}}{6}.$
(52)

Using Eq. 52, we rewrite Eq. 51 to obtain the average deviation of arms away
from the horizontal, a definition similar to $\langle\cos^{2}\theta\rangle$
for fibers, as

$0.50~{}\left(1-\langle p_{3}^{2}\rangle\right)=\frac{1}{30}{S^{t}_{F}}^{-2}$
(53)

Comparing Eqs. 27, 43 , and 53, it is seen that the mean-square deviation of
the orientation from horizontal is the same for fibers, disks and triads when
defined in terms of settling factors based on the inertial rotation of the
respective particles in a quiescent fluid at a 45∘ angle to gravity.

Fig. 3 presents simulation results for the orientational variance of settling
triads obtained by solving the triad velocity, rotation rate and orientational
dynamics using Eqs. 45a and 45b in the Lagrangian stochastic fluid velocity
gradient model. The simulation results are in good agreement with the rapid
settling theory (Eq. 53) for $S^{t}_{F}\gg 1$. Unlike in the case of fibers,
correlations between the transverse orientation of triads and disk-like
particles do not affect the high $S_{F}$ limit, because the transverse
orientation for disk-like particles is very small. Thus, there is no
difference between Lagrangian and rapidly settling particle frame models for
the high $S_{F}$ behavior of disk-like particles.

Fig. 4 compares stochastic simulation results for triads and fibers at
different $Re_{\lambda}$. It may be noted that the general definition of the
settling factor $S_{F}$ based on each particle’s inertial rotation rate nearly
collapses the results for different particle shapes. The results are also
nearly independent of $Re_{\lambda}$.

Figure 3: Mean-squared orientation of triads as a function of settling factor
$S_{F}^{t}$. The squares correspond to simulations, while the lines are
asymptotes derived in the low and high $S_{F}^{t}$ limits. Note that $\theta$
is the angle made by the normal to the triad plane with gravity and hence,
$\langle\cos^{2}\theta\rangle=1$ in the absence of turbulence. The upper
symbols indicate $\langle\cos^{2}\theta\rangle$ and the lower symbols are
$\langle 0.5(1-\cos^{2}\theta)\rangle$, which is the average variance of
deviation of the triad arms away from the horizontal plane.  Figure 4:
Orientational variance of triads and fibers as a function of $S_{F}$ at
different $Re_{\lambda}$. For fibers the orientational variance is
$\langle\cos^{2}\theta\rangle$, while for triads we plot $\langle
0.5(1-\cos^{2}\theta)\rangle$, which is the average variance of deviation of
the triad arms away from the horizontal plane. The transition from isotropic
to preferential alignment happens around the same range of $S_{F}$ as one
might expect. There is a slight difference of the asymptotes, while they share
the $S_{F}^{-2}$ behavior, due to the correlation of the velocity gradient and
orientation vector for fibers in the rapid settling limits in a Lagrangian
stochastic model.

###  2.5 Corrections for Finite Particle Reynolds Number

The theory outlined above can be extended to finite $\textit{Re}_{\ell}$ as
long as $\textit{Re}_{D}\ll 1$ by using the full expressions derived by Khayat
and Cox [31]. However, the resulting non-linear Reynolds number dependency
couples with particle orientation in a non-trivial way when including terms of
order $\mathcal{O}(\ln(\kappa)^{2})$. Solving the force and torque balance
equations is therefore computationally very expensive. Lopez and Guazzelli
[40] suggested a convenient way to handle this case while at the same time
keeping computational cost at a minimum. Since the first-order term in aspect
ratio nicely decouples drag and lift, we can define two constants, $C_{\perp}$
and $C_{R}$, that account for finite Reynolds number and aspect ratio effects:

$\frac{4\pi\mu
LC_{\perp}}{\ln(2\kappa)}\left(\mathbb{1}-C_{R}\boldsymbol{p}\boldsymbol{p}\right)\cdot\boldsymbol{W}-m\boldsymbol{g}=0$
(54)

Here, $C_{\perp}$ accounts for the overall change in drag on a particle
sedimenting at non-zero Reynolds number and $C_{R}$ accounts for the change of
the drag ratio between a particle sedimenting with $\theta=0$ and
$\theta=\pi/2$. In the low Reynolds number limit, $C_{\perp}=1$ and
$C_{R}=1/2$ (see Eq. 3). The expressions for $C_{\perp}$ and $C_{R}$ include
the full analytical expressions given by Khayat and Cox [31] and the only
approximation comes from interpolating the angular dependence at intermediate
orientations. They are defined as:

$C_{\perp}=\frac{\ln(2\kappa)}{\ln(\kappa)}\left(1+\frac{\mathcal{F}_{\perp}}{\ln(\kappa)}\right)$
(55)
$C_{R}=\left(1-\frac{1}{2}\frac{(1-\mathcal{F}_{\perp}/\ln(\kappa))}{(1-\mathcal{F}_{\parallel}/\ln(\kappa))}\right)$
(56)

where $\mathcal{F}_{\parallel}=\mathcal{F}_{D}(\textit{Re}_{\ell},\theta=0)$
and $\mathcal{F}_{\perp}=\mathcal{F}_{D}(\textit{Re}_{\ell},\theta=\pi/2)$.
The above expression has two small parameters present - (ln$(2\kappa))^{-1}$
and (ln$(\kappa))^{-1}$. This is due to the difference in the choice of
perturbation parameters in the slender body theories of [3] and Khayat and Cox
[31]. The expressions for $\mathcal{F}_{D}(\textit{Re}_{\ell},\theta)$ are
given by Khayat and Cox [31]. The two constants $C_{\perp}$ and $C_{R}$ are
plotted as functions of Reynolds number in Fig. 5 for different aspect ratios.
The unexpected behavior of these functions at low Reynolds numbers shows that
the theory is very sensitive to the high aspect ratio requirement and that the
Stokes flow limit can only be recovered when $\kappa\to\infty$ as
$\textit{Re}_{\ell}\to 0$.

Figure 5:  (a) $C_{\perp}$ as a function of Reynolds number
$\textit{Re}_{\ell}$ for three different aspect ratios, $\kappa=20,100$ and
$10^{6}$. The inset shows the $\kappa=20$ results for larger Reynolds number
for comparison with experiments. (b) $C_{R}$ as a function of Reynolds number
$\textit{Re}_{\ell}$ for the same three aspect ratios as in (a). The inset
again shows the $\kappa=20$ results for larger Reynolds number for comparison
with experiments.

Solving Eq. 54 for the velocity of the fiber $\boldsymbol{W}$ yields the same
expression as derived by Lopez and Guazzelli [40], who have approached this
problem in terms of the mobility matrix.

For finite Reynolds numbers, the experimental measurements of the
sedimentation rate of horizontal fibers in quiescent fluid
$W_{\textit{min}}=|\boldsymbol{W}|_{\theta{=}\pi/2}$ described in section 3
enable us to determine $C_{\perp}$. This also includes the uncertainties in
particle dimensions and density. With

$C_{\perp}=\frac{\ln(2\kappa)mg}{4\pi\mu LW_{\textit{min}}},$ (57)

we measure $C_{\perp}\approx 3.5$ and 3.0 for our small fibers and triads, and
$C_{\perp}\approx 7.8$ and 7.6 for our large fibers and triads, respectively.
In order to determine $C_{R}$, one has to measure the velocity of vertical
fibers $W_{\textit{max}}=|\boldsymbol{W}|_{\theta{=}0}$. It is well known that
$W_{\textit{max}}=2~{}W_{\textit{min}}$ for slender fibers in the Stokes flow
limit. Finite aspect ratio effects lower this ratio, but it increases again
slowly when considering finite Reynolds number effects. Our particles fall
into the range where this ratio is approximately 2. As seen in Fig. 5, the
theoretical predictions for $C_{\perp}$ and $C_{R}$ (from Eqs. 55 and 56) are
less than our experimental measurements. This may be due to the finite value
of $\textit{Re}_{D}$ in the experiments. Thus measurements of $C_{\perp}$ and
$C_{R}$ from quiescent fluid experiments are used for the analysis.

Similar to the finite Reynolds number force corrections ($C_{\perp}$ and
$C_{R}$), we introduce a constant $C_{G}$ to correct the inertial torque for
large particles. Here, $C_{G}$ is defined as the ratio of inertial torque from
the full expression $\mathcal{F}_{G}(\textit{Re}_{\ell},\theta)$ (Eq. 6.13)
given by [31] to the low Reynolds number expression
$\mathcal{F}_{G}(\textit{Re}_{\ell}\ll 1,\theta)$ (Eq. 6.22). The behavior of
$C_{G}$ is plotted in Fig. 6 and shows that the inertial torques decrease
quickly with increasing Reynolds number. The value of $C_{G}$ for the
experiments can be determined by comparing the time scales from the low
Reynolds number model Eq. 9 and the empirically determined time scale Eq. 14:

$C_{G}=\frac{\tau_{\textit{sed}}}{T_{\textit{sed}}}$ (58)

and we find $C_{G}=0.007$ for small fibers and triads and $C_{G}=0.002$ for
large fibers and triads. These values are larger than the theoretical value of
$C_{G}$ at the corresponding particle Reynolds numbers and
$\theta=45^{\circ}$, where $C_{G}=0.002$ for small fibers and $C_{G}=0.0002$
for large fibers.

With the force and torque correction factors $C_{\perp},C_{R}$ and $C_{G}$,
the force and torque balance equations for large triads can be written as:

$\sum_{n=1}^{3}\left[-C_{\perp}\left(\mathbb{1}-C_{R}\boldsymbol{p}^{\prime}_{n}\boldsymbol{p}^{\prime}_{n}\right)\cdot\bar{\boldsymbol{W}}^{c}+\hat{\boldsymbol{e}}_{g}\right]=0$
(59a)
$\displaystyle\sum_{n=1}^{3}\left[C_{G}S^{f}_{F}\left(\bar{\boldsymbol{W}}^{c}{\cdot}\boldsymbol{p}^{\prime}_{n}\right)\left(\bar{\boldsymbol{W}}^{c}{\times}\boldsymbol{p}^{\prime}_{n}\right)-4(\mathbb{1}{-}\boldsymbol{p}^{\prime}_{n}\boldsymbol{p}^{\prime}_{n}){\cdot}\bar{\boldsymbol{\Omega}^{\textit{c}}}\right.$
$\displaystyle+\left.4(\boldsymbol{p}^{\prime}_{n}{\times}(\bar{\boldsymbol{\Gamma}}\cdot\boldsymbol{p}^{\prime}_{n}))\right]=0$
(59b)

Figure 6:  $C_{G}$ as function of Reynolds number $\textit{Re}_{\ell}$ for
three different angles, $\theta=0^{\circ},45^{\circ}$ and $90^{\circ}$ (blue,
red and yellow lines, respectively).

##  3 Experiments

In this section, we describe experiments in which we measure the orientation
of particles as they sediment in nearly homogeneous, isotropic turbulence and
compare them with the theoretical predictions from section 2. We also measure
particle motion in quiescent fluid to allow the comparison.

###  3.1 Experimental Methods

Figure 7: a) Autocad (perspective) rendering of the vertical water tunnel.
From bottom to top (color online): inlet chamber for through flow (orange),
pressure plate and honeycomb (green), contraction zone (blue), manifolds
(yellow) for jet array (purple), test section (clear), expansion zone (blue),
flow exit chamber (orange). b) 3D-model of the jet array. c) Side view of a
slice through the jet array shows internal channels, turning vanes and
nozzles.

First we describe the vertical water tunnel which allows control of both the
turbulence intensity and the through flow so that sedimentation can be
balanced with advection to keep particles suspended in the test section for
measurement. Then we describe the imaging and particle fabrication methods.

We identified the non-dimensional parameters, $\textit{Re}_{\lambda}$,
$L/\eta$, $\kappa$, $\rho_{p}/\rho_{f}$, and $S_{F}$ that determine the
sedimentation statistics of non-spherical particles in turbulence. With that
in mind, we constructed a 4.2 m tall, vertical water tunnel (see Fig. 7(a))
that gives us good control over each parameter and allows us to explore a
large range of values. It keeps heavy particles suspended and we can
simultaneously control $S_{F}$ by adjusting the amount of turbulence they
experience. The particles are suspended by a mean through flow, which allows
us to record long, individual trajectories. The density ratio and particle
dimensions, $\rho_{p}/\rho_{f}$, $L$ and $\kappa$, were chosen to yield
particle sedimentation rates that could be supported by the through flow. To
ensure a flat velocity profile of the through flow, the fluid first passes
through a pressure plate (2% open area), a 20 cm tall honeycomb flow
straightener and a contraction zone (area ratio 4:1) before entering the test
section. The exit conditions downstream of the test section are kept symmetric
to the inlet conditions. One of the challenges of the experiments is keeping
the particles suspended without clogging filters, valves or the pump. This is
particularly difficult for fibers, which is the reason why they have been
excluded from the turbulent flow experiments in this paper.

Turbulence is generated and controlled with a jet-array issued from a
3D-printed grid with grid spacing $M=6$ cm and 30% solidity. Using Nylon as
material and selective laser sintering we were able to fabricate the jet array
with internal channels and 40 nozzles (see Fig. 7(b)). Each nozzle can be
triggered independently through a solenoid valve to eject a jet of fluid into
the surrounding flow. In the minimum turbulence configuration, all valves are
closed and the jet-array becomes a passive grid for the through flow. The
streamwise mean fluid velocity $\langle U_{f}\rangle_{z}$ has a small
fluctuating component $u^{\prime}_{z}$, resulting in turbulence intensities
$u^{\prime}_{z}/\langle U_{f}\rangle_{z}$ as small as 7% (see Tab. 1), typical
for passive grid configurations. Here, $u^{\prime}_{(x,y,z)}$, are the
components of the root-mean-square (rms) fluctuating velocity. The system can
also be driven solely through the jet-array, similar to the random jet-array
used by Variano and Cowen [64], achieving much larger turbulence intensities.
The intermediate turbulence regimes can be reached by either adjusting the
number of jets, the duration each jet is firing or the jet velocity. For the
experiments presented here, the turbulence intensity was controlled through
the jet velocity, keeping the average number of jets (8 jets, 20% of the total
number of jets) and the average duration of each jet (1 s $\pm 0.25$ s)
constant. Jets were chosen randomly with some nearest neighbor restrictions.
The generated turbulence is essentially isotropic in the horizontal plane
(span-wise, wall-normal), where $u^{\prime}_{x}\approx u^{\prime}_{y}$, but
has an rms fluctuating velocity in the vertical direction $u^{\prime}_{z}$
(stream-wise) that is about 20% higher than in the horizontal directions as
shown in Tab. 2. Both through flow and jet-array are powered by a 3 hp,
variable speed pump that can produce a through flow of up to 10 cm/s in the
test section and an estimated jet velocity of up to 4 m/s at the nozzles.

| Turb. | Pump Speed | Thru Flow | Total Jet Flow | $u^{\prime}_{z}$ | $\langle U_{f}\rangle_{z}$ | $\langle U_{p}\rangle_{z}$
---|---|---|---|---|---|---|---
| Intensity | [rpm] | [l/s] | [l/s] | [mm/s] | [mm/s] | [mm/s]
Small | 0.07 | 700 | 1.7 | 0 | 1.32 | 19.78 | -2.87
Triads | 0.21 | 800 | 1.3 | 0.4 | 5.01 | 23.38 | 0.22
0.39 | 1200 | 1.4 | 1.0 | 11.64 | 29.98 | 5.48
0.62 | 1550 | 0.9 | 1.5 | 16.68 | 26.91 | 2.26
0.91 | 2000 | 0.2 | 2.1 | 19.98 | 21.96 | -3.98
1.06 | 2100 | 0 | 2.2 | 21.28 | 20.08 | -4.12
Large | 0.07 | 1150 | 2.9 | 0 | 2.39 | 34.14 | -1.92
Triads | 0.10 | 900 | 2.1 | 0.3 | 3.10 | 30.82 | -5.68
0.29 | 1100 | 1.6 | 0.8 | 9.19 | 31.68 | -7.07
0.28 | 1700 | 3.2 | 1.0 | 16.24 | 58.00 | 15.88
0.37 | 1700 | 2.5 | 1.3 | 18.61 | 50.28 | 8.88
0.95 | 2500 | 0.6 | 2.5 | 28.56 | 30.05 | -4.35
Table 1: Experimental parameters: Volumetric through flow and total flow
emitted by the jet array (measured with two magnetic flow meters);
$u^{\prime}_{z}$, stream-wise rms fluctuating velocity;
$\langle\boldsymbol{U}_{f}\rangle$, mean fluid velocity in the detection
volume (measured with tracers). $\langle\boldsymbol{U}_{p}\rangle$, mean
particle velocity.

Two coarse meshes confine the particles in a clear, tall, 30 x 30 x 150 cm3
test section, where four high-speed cameras image a $12$ x $12$ x $10$ cm3
detection volume in the center region, 10$M$ downstream of the jet-array. Two
high-powered, pulsed and monochromatic LED lights (SmartVisionLights ODMOBL
series) create a uniform background illumination. The cameras and the lights
are triggered synchronously at 450 Hz, ensuring a single pulse illumination
during each exposure. A real-time image compression system, which allows
continuous data acquisition for several days, is used to collect large data
sets. It is essential for these experiments, because there is not always a
particle in view on all four cameras, which is required to determine an
accurate particle orientation.

The particles used in the experiments are 3D-printed (using VeroBlack material
and fused deposition modeling, $\rho_{p}=1.13-1.15$ g cm-3) and consist of
several slender fibers, connected at the center (we call them ramified
particles). In this case, three fibers of equal length $L=2\ell$ and radius
$r$ with aspect ratio $\kappa{=}\ell/r{=}20$, oriented in planar symmetry and
with a 120∘ angle between them, form a triad (see Fig. 1). The triad is a
model for a very small aspect ratio disk-like particle. The advantage of using
ramified particles as models is that the orientations can be measured very
accurately and even rotations around the symmetry axis can be resolved. The
rotations of any ellipsoid can be approximated by a corresponding ramified
particle (exact for small particles) and therefore even the rotations of
spherical particles can be tracked accurately, overcoming one measurement
limitation inherent to ellipsoidal particles.

Small Triads
---
Turb. | $\textit{Re}_{\lambda}$ | $u^{\prime}_{(x,y,z)}/u^{\prime}_{z}$ | $\bar{u}$ | $\mathcal{L}$ | $\epsilon$ | $\eta$ | $\tau_{\eta}$
Intensity |  |  | [mm/s] | [mm] | [mm2/s3] | [mm] | [s]
0.07 | 29 | (0.81, 0.83, 1.00) | 1.38 | 45 | 0.04 | 2.16 | 5.11
0.21 | 91 | (0.83, 0.80, 1.00) | 4.41 | 114 | 0.75 | 1.00 | 1.12
0.39 | 141 | (0.85, 0.85, 1.00) | 10.34 | 116 | 9.5 | 0.53 | 0.31
0.62 | 162 | (0.83, 0.83, 1.00) | 14.80 | 108 | 30 | 0.39 | 0.18
0.91 | 192 | (0.87, 0.87, 1.00) | 18.31 | 123 | 50 | 0.35 | 0.14
1.06 | 194 | (0.85, 0.86, 1.00) | 19.25 | 119 | 60 | 0.33 | 0.13
Large Triads
0.07 | 36 | (0.80, 0.83, 1.00) | 2.11 | 38 | 0.25 | 1.32 | 1.91
0.10 | 56 | (0.82, 0.80, 1.00) | 2.74 | 69 | 0.30 | 1.26 | 1.77
0.29 | 102 | (0.86, 0.86, 1.00) | 8.32 | 77 | 7.5 | 0.56 | 0.35
0.28 | 153 | (0.85, 0.85, 1.00) | 14.39 | 99 | 30 | 0.40 | 0.17
0.37 | 162 | (0.86, 0.84, 1.00) | 16.83 | 95 | 50 | 0.35 | 0.14
0.95 | 200 | (0.86, 0.85, 1.00) | 25.75 | 95 | 180 | 0.25 | 0.07
Table 2: Turbulence parameters.
$\textit{Re}_{\lambda}=\sqrt{15\bar{u}\mathcal{L}/\nu}$, Reynolds number;
$u^{\prime}_{(x,y,z)}$, components of the rms fluctuating velocity normalized
by $u^{\prime}_{z}$; $\bar{u}=\langle u^{\prime}_{i}\rangle_{i}$, component-
average rms fluctuating velocity; $\mathcal{L}=\bar{u}^{3}/\epsilon$, energy
input length scale; $\epsilon$, mean energy dissipation rate;
$\eta=(\nu^{3}/\epsilon)^{1/4}$, Kolmogorov length;
$\tau_{\eta}=(\nu/\epsilon)^{1/2}$, Kolmogorov time.

We fabricated 150 triads with smallest dimension $r=(225\pm 5)$ $\mu$m
(referred to as small particles, but not small enough to be in the Stokes flow
regime) and with $r=(450\pm 5)$ $\mu$m (large particles), both with
$\kappa=20$. The particle Reynolds number $\textit{Re}_{D}=wD/\nu$ based on
the fiber diameter ranges from $\textit{Re}_{D}\approx 10$ for small triads to
$\textit{Re}_{D}\approx 40$ for large triads (see Table 3), where $w$ is the
velocity of the particle relative to the fluid. The fluid viscosity is
$\nu=0.9131\pm 0.02\times 10^{-6}$ m2 s-1, with the uncertainty coming from
temperature fluctuations ($T=24.0^{\circ}\pm 1^{\circ}$). Fibers used for some
parts of the experiments were manually cut from Nylon fishing line with a very
similar density of $\rho_{p}=1.13-1.15$ g cm-3, but much smoother surface than
the 3D-printed triads. The fiber radius and aspect ratio was close to that of
the arms that make up the triads. Grey, neutrally buoyant micro spheres with a
diameter of 250 $\mu$m were used as tracer particles to measure the fluid
velocity, calculate structure functions and extract mean energy dissipation
rates $\epsilon$. For the experiments with small triads, the local fluid
velocity at the particle position was measured using tracers within a sphere
of radius $3L$ around the particle. The local mean fluid velocity
$\langle\boldsymbol{u}_{f}\rangle$ did not depend very strongly on the size of
the tracer-sphere and a relatively large radius was chosen to minimize
measurement noise. This was not successful for the experiments with the large
triads due to an insufficient number of tracers and so the overall mean fluid
velocity $\langle\boldsymbol{U}_{f}\rangle$ was used to calculate a relative
particle velocity. In our experiments, the non-dimensional concentration is
$nL^{3}=O(10^{-3})$, where $n$ is the number density, and thus the effect of
hydrodynamic interactions is negligible.

The particles have to be suspended near the center of the test section in
order to take continuous data and gather enough statistics. Depending on
particle size, the through flow can be adjusted to match the particles
sedimentation rate in quiescent fluid. As we increase the turbulence
intensity, triads start to rotate around their equilibrium sedimentation
orientation which, in return, increases their sedimentation rate. The through
flow was adjusted to keep as many triads suspended as possible. For the
largest turbulence intensities, triads were almost entirely suspended by
strong jets from the jet array, whereas for intermediate turbulence
intensities, the through flow had to be increased to lift particles up into
the detection volume. Triads exposed to the lowest turbulence intensities show
the same orientation and sedimentation statistics as triads in quiescent
fluid. The volumetric mean flow and the total volumetric flow emitted by the
jet array were measure with two separate magnetic flow meters (Toshiba GF630
series), see Table 1.

###  3.2 Sedimentation in Quiescent Fluid

Before presenting experimental results on the orientation of particles
sedimentation in turbulent flows, we will first document the translational and
orientational motion of the particles in quiescent fluid. The Reynolds numbers
based on both rod diameter and arm length are larger than the range for which
the theoretical calculations were performed. Thus, we will use these
measurements to quantify the rate of rotation of the triads in the experiments
and define the settling parameter $S_{F}$. The experiments will also allow an
assessment of the extent to which the orientation dependence of the settling
velocity and rotation rate resembles that of low Reynolds number particles so
that a theory based on adjusted values of $S_{F}$ might capture the
orientation of particles in turbulent flows.

The experiments in quiescent fluid include two kinds of particles, fibers and
triads, with two different sizes. Fibers were chosen to match the non-
dimensional parameters of triads as closely as possible. In quiescent fluid,
both particles will assume a stable sedimentation orientation with its longest
axis perpendicular to its sedimentation direction. In the lab reference frame,
$\boldsymbol{\hat{z}}$ is upwards and gravity $\boldsymbol{g}$ is downwards,
this means the stable orientation of fibers is $p_{z}=0$ ($\theta=90^{\circ}$)
and $p_{z}=1$ ($\theta=0$) for triads ($p_{z}=p_{3}$).

The measurements of $W_{\textit{min}}$ together with Eq. 57 determine
$C_{\perp}$ of the particles in their stable orientation. Moreover, the
orientation distributions and variances of particles, sedimenting in quiescent
fluid, contain valuable information about particle inhomogeneities and
fabrication defects.

To gain insight into the sedimentation statistics in quiescent fluid, at
orientations other than their equilibrium orientation, we disturb the
particles by letting them hit a thin nylon string and recording the resulting
trajectories. The sedimentation rate is measured for all orientations along
the particles trajectory and therefore includes effects of the particle’s
history. In other words, the sedimentation rate is measured during a transient
phase whereas the ramified particle model assumes quasi-steady state
sedimentation velocity.

Figure 8: a) Components of the relative particle velocity $\boldsymbol{W}$ as
function of particle orientation, $p_{z}$, measured in quiescent fluid. The
top curves show the vertical component $W_{z}$ for small fibers and triads
(small and large circles and triangles, respectively). The lower curves show
the horizontal component $W_{h}$, where
$\boldsymbol{\hat{h}}=\boldsymbol{\hat{z}}\times(\boldsymbol{p}\times\boldsymbol{\hat{z}})$.
Both are normalized by the sedimentation velocity of a broadside settling
particle, $W_{\textit{min}}$. Error bars are showing the standard deviation
between individual trajectories. The dashed lines show the predictions of the
simple ramified particle model for infinitely long fibers (thin disks) with
$M_{\perp}{=}2M_{\parallel}$. b) The long axis of a fiber with $p_{z}{=}0.5$
is oriented at $30^{\circ}$ with respect to the horizontal. The inset shows a
photograph of the smooth surface of the nylon fibers. c) The plane of a triad
with $p_{z}{=}0.5$ is oriented at $60^{\circ}$ with respect to the horizontal.
The inset shows a photograph of a 3D printed triad with a rougher surface and
features up to 0.1$r$, for $r=0.225$ mm.

In their equilibrium orientation, our small fibers and triads sediment at a
mean rate of $W_{\textit{min}}=19.8$ mm/s and 23.2 mm/s, respectively, and our
large fibers and triads at $W_{\textit{min}}=35.9$ mm/s and 36.8 mm/s,
respectively.

Figure 8 (a) shows the components of the relative particle velocity
$\boldsymbol{W}$ of fibers and triads as functions of particle orientation.
The top curves show the mean vertical component $W_{z}$. Fibers of both sizes
follow the predictions of the theoretical model from section 3.4 very well,
even though they are far outside the range of Reynolds numbers where this
model is valid. Their sedimentation rate in the vertical orientation is about
twice their sedimentation rate in the horizontal orientation, $\langle
W_{z}|_{p_{z}{=}1}\rangle\approx 2~{}\langle W_{z}|_{p_{z}{=}0}\rangle$.
Triads on the other hand are not quite reaching the ratio of sedimentation
rates (1.5) predicted by the ramified particle model, but $\langle
W_{z}|_{p_{z}{=}0}\rangle\approx 1.4~{}\langle W_{z}|_{p_{z}{=}1}\rangle$.
There are multiple reasons that could explain this lower ratio of
sedimentation rates. For our particles, the low fiber-diameter Reynolds number
and high aspect ratio approximations are not valid ($\textit{Re}_{D}>10$ and
$\kappa=20$) and both effects are known to change that ratio [31]. As
described before, one can model these effects in the simple ramified particle
model by adjusting the constant $C_{R}$. However, if this was the dominant
reason for a lower ratio of sedimentation rates, we would expect to see a
similar effect for fibers. The most likely reason why we see a discrepancy
between the model prediction and the measurements of triads is that the model
neglects the interactions between arms of a ramified particle. We also refer
the interested reader to Chapter 4 of Kramel [36] for a comparison of the
angular dependence of the relative velocity of the triad and fluid obtained
from theory and experiments.

In addition to vertical sedimentation, fibers and triads have a non-zero
horizontal settling velocity depending on particle orientation. The bottom
curves in figure 8 (a) show the mean horizontal component,
$W_{h}=\boldsymbol{W}\cdot\boldsymbol{\hat{h}}$, where
$\boldsymbol{h}=(\mathbb{1}-\boldsymbol{\hat{g}}\boldsymbol{\hat{g}})\cdot\boldsymbol{p}$
is the projection of $\boldsymbol{p}$ in the plane perpendicular to
$\boldsymbol{g}$ (horizontal). Based on the model, we expect this component to
be maximized independently of shape when $\theta=45^{\circ}$, or
$p_{z}=1/\sqrt{2}$. The experiments show that both fibers and triads reach the
largest horizontal velocity at a different angle, when $p_{z}\approx 0.5$. One
has to keep in mind that $\boldsymbol{p}$ points along the symmetry axis of
the fiber, but is perpendicular to the plane of the triad and therefore
$p_{z}=0.5$ for fibers means the long axis makes an angle of $30^{\circ}$ with
respect to the horizontal, while for triads the long axis makes an angle of
$60^{\circ}$ with respect to the horizontal, see Fig. 8 (b) and (c). The
experiments show that fibers have roughly twice the horizontal velocity of
triads for all orientations, which is in agreement with the ramified particle
model.

Figure 9: Angle between the force of gravity and particle orientation as
function of time. Fibers (blue circles) approach their stable orientation
where $\boldsymbol{p}$ is perpendicular to $\boldsymbol{\hat{g}}$, whereas
triads (red triangles) are stable when $\boldsymbol{p}$ is parallel to
$\boldsymbol{\hat{g}}$. The black circles show the theoretical prediction for
small fibers from the simulation results of Shin et al. [58]. a) Normalized by
the inverse of the predicted rotation-rate at 45 degree using the tumbling
rate due to inertial torque (Eq. 8)

. b) Normalized by the inverse of the measured rotation-rate at 45 degree.

The model for inertial torques, Eq. 4, is valid at low particle Reynolds
numbers. Shin et al. [58] have shown in simulations, that at
$\textit{Re}_{\ell}\sim 10$, the inertial torque passes through a maximum and
begins decreasing with increasing Reynolds number. We can therefore expect
that the time scale $\tau_{\textit{sed,}45}$, predicted by the low Reynolds
number model (Eq. 9) is too short for the particles in the experiments. Figure
9 (a) shows the angle between gravity $\boldsymbol{g}$ and the particle
orientation $\boldsymbol{p}$ as function of time, normalized by
$\tau_{\textit{sed,}45}$. The trajectories shown are averaged trajectories of
many individual fibers and triads, with $t=0$ chosen when each particle is at
$\theta=45^{\circ}$. The error bars are showing the standard deviation between
individual trajectories, which is zero at $t=0$ because we choose
$\theta=45^{\circ}$ to define time zero. We also show the simulation data from
Shin et al. [58] of fibers with $Re_{D}\ll 1$. For our particles, the low
Reynolds number model clearly over-estimates the strength of the inertial
torques. Furthermore, the low Reynolds number time scale for particle
rotation, $\tau_{\textit{sed,}45}$, does not collapse the experimental curves
with one another nor the experiments with the simulations. For that reason, we
can not use the definition of $\tau_{\textit{sed,}45}$ from the low Reynolds
number model to calculate $S_{F}$ for the experimental particles with
$\textit{Re}_{D}=10$.

To gain information about the strength of the inertial torques when the
particle Reynolds number is no longer small, we extract the time it takes a
particle to come to alignment with its equilibrium orientation. We measure the
the rotation rate of the particles when $\theta=45^{\circ}$ and use it to
define an empirical time scale of the inertial torques as
$T_{\textit{sed,}45}=1/|\boldsymbol{\dot{p}}|_{\theta=45^{\circ}}$. Here, we
determine the rotation rate when $\theta=45^{\circ}$ by fitting a straight
line to the measurements of $p_{z}$ over the range $t=0$ to
$0.05T_{\textit{sed,}45}$. Figure 9 (b) uses this definition to collapse the
curves. It is notable that the dependence of angle on $t/T_{\textit{sed}}$ is
similar in the experiments and theory suggesting that using the measured
$T_{\textit{sed}}$ the theory may predict the particle dynamics well.
Interestingly, $T_{\textit{sed,}45}$ is roughly the same for all the particle
sizes and shapes used in our experiments with $T_{\textit{sed,}45}=(1.8\pm
0.1)$ s for small fibers and $T_{\textit{sed,}45}=(1.9\pm 0.1)$ s for large
fibers and $T_{\textit{sed,}45}=(1.7\pm 0.1)$ s for small triads and
$T_{\textit{sed,}45}=(1.9\pm 0.1)$ s for large triads. One notable difference
between the predicted and observed orientational dynamics is that triads often
significantly overshoot the $\theta=0$ point (not shown) and return at
slightly different rates. Averaging this over different trajectories makes it
appear that the equilibrium angle is at $\theta\sim 10^{\circ}$ in the shown
time range. In the long time limit $\theta$ will approach $0^{\circ}$. We will
use this time scale of the inertial torques, $T_{\textit{sed,}45}$, to
calculate an empirical value of the settling factor $S_{F}$ for the turbulence
experiments.

$\kappa$ | $\rho_{p}/\rho_{f}$ | $\textit{Re}_{\lambda}$ | $L/\eta$ | $S_{F}$ | $\textit{Re}_{D}$
---|---|---|---|---|---
20 | 1.146 | 29 (36) | 4.5 (14.3) | 1.95 (1.30) | 10.82 (35.74)
91 (56) | 9.0 (14.3) | 0.66 (1.17) | 11.06 (35.26)
141 (102) | 17.0 (32.1) | 0.26 (0.44) | 11.88 (36.89)
162 (153) | 23.1 (45.0) | 0.18 (0.25) | 12.31 (40.34)
192 (162) | 25.7 (51.4) | 0.15 (0.22) | 12.36 (39.58)
194 (200) | 27.3 (72.0) | 0.14 (0.15) | 12.21 (33.25)
Table 3: Non-dimensional parameters of small (and large) triads. $\kappa=L/D$,
aspect ratio; $\rho_{p}/\rho_{f}$, density ratio;
$\textit{Re}_{\lambda}{=}\sqrt{15\bar{u}\mathcal{L}/\nu}$, Reynolds number;
$L/\eta$, non-dimensional particle size, with length $L=9$ mm ($L=18$ mm);
$S_{F}$, settling number; . $\textit{Re}_{D}$, particle Reynolds number based
on the diameter .

###  3.3 Sedimentation in Turbulence

The sedimentation of ramified particles under different turbulence intensities
gives insight into the competition between turbulence, which tend to randomize
particle orientations, and inertial torques, which align particles with their
stable sedimentation orientation. We present orientation distributions of
triads as the variance of $p_{z}$ for various values of the turbulence
intensity and particle size. We also compare the experiments to simulations
and theoretical predictions based on slender body theory.

The relative particle orientation can be defined using two angles, $\theta$
and $\psi$. We define these two angles as

$\displaystyle\cos(\theta)$
$\displaystyle=\left|\boldsymbol{p}\cdot\boldsymbol{\hat{g}}\right|$ (60)
$\displaystyle\cos(\psi)$
$\displaystyle=\left|\frac{\left[\left(\mathbb{1}-\boldsymbol{p}\boldsymbol{p}\right)\cdot\boldsymbol{\hat{g}}\right]\cdot\boldsymbol{p^{\prime}}}{|\left(\mathbb{1}-\boldsymbol{p}\boldsymbol{p}\right)\cdot\boldsymbol{\hat{g}}|}\right|$
(61)

where $\theta$ is the angle between $\boldsymbol{p}$ and gravity and $\psi$
the angle between an arm $\boldsymbol{p}^{\prime}$ and the projection of
gravity into the plane of the particle
$\left(\mathbb{1}-\boldsymbol{p}\boldsymbol{p}\right)\cdot\boldsymbol{\hat{g}}$.
In isotropic turbulence, the third Euler angle is a rotation around
$\boldsymbol{\hat{z}}$ which is randomly distributed and does not encode any
relevant statistics.

The orientation of an arbitrary rigid body would require the specification of
the three Euler angles, $\theta,\phi,\psi$, and thus we would define a PDF
$\Pi(\theta,\phi,\psi)$. For an isotropic flow, the PDF would be independent
of $\phi$ $\Pi(\theta,\phi,\psi)=\mathcal{P}(\theta,\psi)/2\pi$. Further for
fore-aft symmetric particles a measure of the distribution of the normal to
the plane of the particle can be defined as $P(p_{z})$ -

$\displaystyle
P(p_{z}=\cos\theta)=2\int_{0}^{2\pi}\mathcal{P}(\theta,\psi)d\psi.$ (62)

Thus as per our definition, for isotropic distribution, $P(p_{z})=1$. The
current theoretical approach does not depend on $\psi$, but the experimental
results display $\psi$ dependence. This behavior is likely due to inertial
effects at larger Reynolds numbers that are not included in the theory and
possibly some gravitational torques due to differences in the arms.To observe
the $\psi$ dependence in experiments, we define a PDF

$\displaystyle\Psi(\psi)=\int_{0}^{\pi}\mathcal{P}(\theta,\psi)\sin\theta
d\theta.$ (63)

Figure 10 shows how $P(p_{z})$ changes with turbulence intensity. For low
turbulence intensities, small and large triads are strongly aligned with the
direction of gravity, within a few degrees, reflected in the sharp peak near
$p_{z}=1$. This alignment becomes weaker as the turbulence intensity
increases. The orientation PDFs become more uniform. Even for high turbulence
intensities, particle orientations are not fully randomized. The corresponding
settling factors are summarized in Tab. 4.

Figure 10:  PDF $P(p_{z})$ of the particle orientation at different turbulence
intensities from low to high (colormap cold to hot). $|p_{z}|=1$ means
horizontal alignment, $|p_{z}|=0$ vertical alignment and the PDF of randomly
oriented particles is uniform at 1. (a) Experiments with small triads. (b)
Experiments with large triads. Figure 11:  PDF $\Psi(\psi)$ of the particle
orientation within the plane spanned by the arms of the particle at different
turbulence intensities from low to high (colormap cold to hot). $\psi=0$ means
an arm of the particle is aligned parallel with the direction of gravity
($120^{\circ}$ symmetry), $\psi=\pi/3$ mean an arm of the particle is anti-
parallel to the direction of gravity. (a) Experiments with small triads. (b)
Experiments with large triads.

Triads also show preferential alignment within the plane of the particle.
Figure 11 shows the PDF $\Psi(\psi)$. It is surprising that for low turbulence
intensities, small triads show a strong peak at $\psi=\pi/3$, meaning any one
of the three arms is pointing slightly upward. We assume that particle defects
or fabrication inhomogeneities cause this alignment, e.g. one arm might
experience slightly larger drag. We do not see such behavior for large triads,
where these effects would have a smaller impact. With increasing turbulence
intensity, the particles get kicked out of their equilibrium orientation and
interactions between the arms cause one arm to preferentially align parallel
to the direction of gravity, pointing downward. Large triads show that this
effect is strongest for intermediate turbulence intensities, where turbulent
fluctuations are strong enough to kick a particle far enough out of its
equilibrium orientation so that interactions become relevant, but do not fully
randomize the particles orientation yet. The theoretical model which assumes
symmetric arms and neglects hydrodynamic interactions among the arms predicts
a uniform distribution of $\psi$.

Figure 12 shows spherical histograms of particle orientations for small and
large triads for different values of $S_{F}$. The histogram depicts the
components of a unit vector defined in spherical coordinates by setting the
polar angle equal to $\theta$, and the azimuthal angle equal to $\psi$. The
large probability at the poles (see Fig. 12 (a) small triads and Fig. 12 (d)
large triads at the lowest turbulence intensity) indicates that
$\boldsymbol{p}$ is strongly aligned with the direction of gravity. The
symmetric, but not random probability distribution around the equator (120∘
symmetry) shows the preferential alignment of one arm parallel to the
direction of gravity. This preferential alignment can only be seen when the
particles are significantly kicked out of their equilibrium orientation, as
mentioned before. In fact, for small triads at the lowest turbulence
intensity, the probability distribution shows small peaks, indicating opposite
alignment with one arm upward.

Figure 12: Orientation Probability Distribution Function
$\mathcal{P}(\theta,\psi)$ for small triads at (a) $S_{F}=1.95$, (b)
$S_{F}=0.18$, and (c) $S_{F}=0.14$ and large triads at (d) $S_{F}=1.30$, (e)
$S_{F}=0.25$, and (f) $S_{F}=0.15$.  Figure 13:  Mean-square particle
orientation as a function of the settling factor $S_{F}$. Experimental results
for large and small triads are shown as large and small diamonds with dashed
lines, respectively, while the simulations are circles. The horizontal line
indicates isotropic orientation. Note that the definition of $S_{F}$ for the
experimental results uses the experimentally measured inertial rotation rate
of a triad at a 45∘ angle to gravity in a quiescent fluid and the eddy
turnover time for eddies of the particle size.

To characterize the alignment of the normal vector to the plane of the
particle relative to gravity, we present the variance of the angle of the
triad arms, $\langle p_{z}^{\prime
2}\rangle=0.5\left(1-\langle\cos^{2}(\theta)\rangle\right)$, as a function of
the settling factor $S_{F}$ in Figure 13. It can be seen that simulations and
the experiments with both triad sizes exhibit transitions from a nearly
isotropic distribution corresponding to $\langle p_{z}^{\prime
2}\rangle=0.5\left(1-\langle\cos^{2}(\theta)\rangle\right)=1/3$ at small
$S_{F}$ to a nearly horizontal orientation $\langle p_{z}^{\prime
2}\rangle=0.5\left(1-\langle\cos^{2}(\theta)\rangle\right)\rightarrow 0$ at
large $S_{F}$. The transition occurs over approximately the same range
$S_{F}=0.2$ to $2$ of settling factors for the two particles sizes and for the
theory. This suggests that accounting for the measured rotation rate of the
particles in quiescent fluid and using a filtered turbulent velocity gradient
based on the particle size in defining $S_{F}$ captures the gross features of
particle alignment successfully. At intermediate $S_{F}$ values between about
0.1 and 0.5, the experimental results for the deviation from perfect alignment
are generally somewhat lower than the simulation results. One possible cause
of this difference is that the simulations use a stochastic turbulent velocity
gradient in a Lagrangian reference frame whereas the velocity gradients seen
by particles are decorrelated by particle settling as well as turbulence
evolution. For $S_{F}\gg 1$, the rapid settling theory predicts that $\langle
p_{z}^{\prime 2}\rangle=0.5\left(1-\langle\cos^{2}(\theta)\rangle\right)$ is
proportional to $S_{F}^{-2}$ and the simulations follow this scaling. The
steep decline of the experimental orientation variance from $S_{F}=0.5$ to 2
is consistent with evolution toward this limiting behavior. However, the
velocity variance for the highest $S_{F}$ measurement for the small triads is
clearly trending above this limit. This is likely due to imperfections in the
particle fabrication resulting in a slight deviation from horizontal alignment
that was seen even in the quiescent fluid experiments.

Small Triads |  | Large Triads
---|---|---
$S_{F}$ | $\langle\cos^{2}(\theta)\rangle$ | $\frac{1}{2}(1-\langle\cos^{2}(\theta)\rangle)$ |  | $S_{F}$ | $\langle\cos^{2}(\theta)\rangle$ | $\frac{1}{2}(1-\langle\cos^{2}(\theta)\rangle)$
1.95 | 0.970 | 0.015 |  | 1.30 | 0.982 | 0.009
0.66 | 0.933 | 0.034 |  | 1.17 | 0.981 | 0.010
0.26 | 0.757 | 0.122 |  | 0.44 | 0.853 | 0.074
0.18 | 0.606 | 0.197 |  | 0.25 | 0.559 | 0.221
0.15 | 0.569 | 0.216 |  | 0.22 | 0.534 | 0.233
0.14 | 0.554 | 0.223 |  | 0.15 | 0.431 | 0.285
Table 4: Sedimentation parameters. Settling number defined empirically for
triads, $S_{F}$; $\langle\cos^{2}(\theta)\rangle$ and
$0.5(1-\langle\cos^{2}(\theta)\rangle)$, orientation variance of
$\boldsymbol{p}$ and one of the arms $\boldsymbol{p}^{\prime}$.

####  3.3.1 A model for orientation PDF - non-Gaussian effects

The variation of the orientation PDF with the turbulence intensity (figure 10)
highlights the non-Gaussian nature of the distributions. Anand et al. [1]
calculated the higher moments of the orientation PDFs from their DNS
calculations and showed that the orientation PDFs are non-Gaussian due to the
non-Gaussian nature of the turbulent velocity gradient. Here we consider a
model problem to analyze the non-Gaussian orientation PDFs of anisotropic
particles settling in a turbulent flow. In our model problem, we will consider
the rapid settling of thin disks in vertical simple shear flows (flow axis
aligned with gravity) as shown in figure 14. We have seen earlier that the
orientation dynamics of disks closely resemble that of triads. Our objective
here is to obtain the dependence of the equilibrium orientation of the disks
on the shear rate, which would then help us derive the orientation PDF of the
disks based on the assumed form of PDF for the shear rates.

Figure 14: A thin disk settling in a vertical shear flow

For a vertical simple shear flow ($\Gamma_{ij}=\gamma\delta_{3i}\delta_{1j}$),
the evolution equation for orientation (see equation 37) reduces to

$\displaystyle\dot{\theta}=-\gamma\cos^{2}\theta-\frac{3cW_{3}^{2}}{2\nu}\sin\theta\cos\theta$
(64)

The above system has two fixed points - an unstable one at $\theta=\pi/2$ and
a stable one at $\theta=-\tan^{-1}\left(2\gamma\nu/(3cW_{3}^{2})\right)$. Thus
at equilibrium, we have $\tan\theta=-1/S_{F}^{\textrm{SS}}$. Motivated by our
earlier calculations, we have introduced a settling parameter for a thin disk
settling in a vertical simple shear flow
$S_{F}^{\textrm{SS}}=2\gamma\nu/(3cW_{3}^{2})$. The relationship between the
equilibrium angle and the shear rate allows us to construct the PDF for
$\theta$ or $p_{z}$ with a known PDF for $\gamma$.
In our earlier calculation of the Lagrangian model for velocity gradient, we
used the Girimaji and Pope [19] model. Their model incorporates the lognormal
distribution for pseudodissipation. To construct a PDF for the shear rate,
$P_{\gamma}$, we propose that $\gamma=h\bar{\phi}^{1/2}$, where $h$ is a
normalized shear rate that obeys a Gaussian distribution with unit variance
($P_{h}(h)$) and $\bar{\phi}$ is proportional to the pseudodissipation and is
lognormally distributed. The variance of the logarithm of $\bar{\phi}$,
$\sigma^{2}_{\log\bar{\phi}}$, depends on $Re_{\lambda}$ and is obtained from
expressions provided by Koch and Pope [35]. Using the known PDFs, $P_{h}$ and
$P_{\bar{\phi}}$, we obtain

$\displaystyle
P_{\gamma}\left(|\gamma|\right)=4k|\gamma|\exp\left[-\frac{2(\log|\gamma|)^{2}}{\sigma_{\log\bar{\phi}}}\right]\int_{0}^{\infty}P_{h}(h)\exp\left[-\frac{2(\log|h|)^{2}}{\sigma_{\log\bar{\phi}}}\right]h^{\frac{4\log|\gamma|}{\sigma_{\log\bar{\phi}}}-2}dh,$

where $k=e^{-\sigma_{\log\bar{\phi}}/2}/\sqrt{2\pi\sigma_{\log\bar{\phi}}}$.
We next use the relation between the equilibrium angle and shear rate to
obtain the PDF $P(p_{z})$ from $P_{\gamma}$. To compare this model for the
orientation distribution with experiments, we choose a value
$S_{F}^{\textrm{SS}}$ of the settling parameter that gives the same
orientational variance as the experimental measurements

$\displaystyle\langle\cos^{2}\theta\rangle_{\textrm{model}}=\langle\cos^{2}\theta\rangle_{\textrm{expts.}}.$
(66)

Figure 15: Comparison of PDFs $P(p_{z})$ between experiments and model for
small triads. Red symbols and line correspond to $S_{F}=0.66$ and
$S_{F}^{\textrm{SS}}=6.76$ respectively. Blue symbols and line correspond to
$S_{F}=0.26$ and $S_{F}^{\textrm{SS}}=2.59$ respectively. $Re_{\lambda}$ is
the same in both experiments and model - $Re_{\lambda}=(91,141)$. Comparison
of higher moment - From the experiments
$\langle(1-p_{3}^{2})^{2}\rangle\approx(0.0119,0.1593)$ and from the model
$\langle(1-p_{3}^{2})^{2}\rangle=(0.0159,0.1215)$

.

In figure 15, we show comparisons of PDF for $p_{z}$ from our model problem
with that obtained from experiments. The orientation PDFs obtained from the
model problem compare well with those from experiments, with the comparisons
getting better for higher $S_{F}^{\textrm{SS}}$ scenarios. We also compare the
higher moments, $\langle(1-p_{3}^{2})^{2}\rangle$, and as visible in figure
15, they agree reasonably with those computed from the experimental data. The
model also allows exploration of the dependence of the orientation PDF on
$Re_{\lambda}$, which modulates the non-Gaussian nature of the velocity
gradient statistics.

(a)

(b)

Figure 16: Variation of the orientation moments calculated from the model
$P(p_{z})$ compared with their large $S_{F}^{\textrm{SS}}$ asymptotic forms
(equation 71). The inset of figure (b) shows the variation of the kurtosis
with $Re_{\lambda}$ in the $S_{F}^{\textrm{SS}}\gg 1$ limit.

The rapid settling analysis outlined above will also allow for the evaluation
of particle orientation tensors, $p_{i}p_{j}$ and $p_{i}p_{j}p_{k}p_{l}$, that
would be required to find the stresslet. Although the particle orientation is
highly correlated with the local flow, the orientational moments obtained here
are averaged over all flows. The particle orientation tensors can be written
as,

$\displaystyle\langle p_{i}p_{j}\rangle$ $\displaystyle=$
$\displaystyle\delta_{ij}+\lambda\left(\delta_{i3}\delta_{j3}-\frac{1}{3}\delta_{ij}\right),$
(67) $\displaystyle\langle p_{i}p_{j}p_{k}p_{l}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{5}\left(1-\bar{\lambda}-\frac{\lambda}{3}\right)\left(\delta_{ij}\delta_{kl}+\delta_{ik}\delta_{jl}+\delta_{il}\delta_{jk}\right)+\left(\lambda-7\bar{\lambda}\right)\delta_{i3}\delta_{j3}\delta_{k3}\delta_{l3}+$
$\displaystyle\bar{\lambda}\left(\delta_{i3}\delta_{j3}\delta_{kl}+\delta_{i3}\delta_{k3}\delta_{jl}+\delta_{i3}\delta_{l3}\delta_{jk}+\delta_{j3}\delta_{k3}\delta_{il}+\delta_{j3}\delta_{l3}\delta_{ik}+\delta_{k3}\delta_{l3}\delta_{ij}\right),$

where,

$\displaystyle\lambda$ $\displaystyle=$ $\displaystyle\frac{3}{2}\left(\langle
p_{3}^{2}\rangle-1\right)=\langle\textrm{P}_{2}\left(p_{3}\right)\rangle-1,$
(69) $\displaystyle\bar{\lambda}$ $\displaystyle=$
$\displaystyle\frac{6\langle p_{3}^{2}\rangle-5\langle
p_{3}^{4}\rangle-3}{8}=\frac{\langle\textrm{P}_{2}\left(p_{3}\right)\rangle-\langle\textrm{P}_{4}\left(p_{3}\right)\rangle}{7}-\frac{1}{4},$
(70)

where $\textrm{P}_{n}\left(p_{3}\right)$ is the Legendre polynomial of order
$n$. Rapid settling theory developed in the current study provides us with the
expression for $\langle p_{3}^{2}\rangle$ (equation 42). $\langle
p_{3}^{4}\rangle$ requires information regarding the fourth moment of the
turbulent velocity gradient tensor. The fourth moment,
$\langle\Gamma_{ij}\Gamma_{kl}\Gamma_{mn}\Gamma_{pq}\rangle$, is a measure of
the non-Gaussian statistics. It involves 105 isotropic tensors constrained by
four invariants [62] that are obtained from DNS or experiments (see Fang et
al. [17] for further details). $\langle p_{3}^{2}\rangle$ and $\langle
p_{3}^{4}\rangle$ can be obtained from the experiments. Table 4 lists the
values of $\langle p_{3}^{2}\rangle$ from the current experiments. The model
$P(p_{z})$ developed in this section (from equation 3.3.1) also allows us to
calculate the values of the orientation moments and the scalar constants,
$\lambda$ and $\bar{\lambda}$, for the particle orientation tensor. Besides
evaluating the moments for arbitrary values of $S_{F}^{\textrm{SS}}$ we can
also find the asymptotic forms in the rapid settling limit,

$\displaystyle\langle(1-p_{3}^{2})^{n}\rangle\sim\frac{\mathcal{C}_{n}(Re_{\lambda})}{(S_{F}^{\textrm{SS}})^{2n}},$
(71)

where
$\mathcal{C}_{n}(Re_{\lambda})=\int_{0}^{\infty}\gamma^{2n}P_{\gamma}\left(\gamma\right)d\gamma$
is a constant that can be found numerically. Figure 16 shows the variation of
these statistical quantities with the settling parameter
($S_{F}^{\textrm{SS}}$) for two values of $Re_{\lambda}$, compared with the
large $S_{F}^{\textrm{SS}}$ asymptotes. In the $S_{F}^{\textrm{SS}}\gg 1$
limit, we can also calculate the kurtosis

$\displaystyle\lim_{S_{F}^{\textrm{SS}}\rightarrow\infty}\frac{\langle(1-p_{3}^{2})^{4}\rangle}{\langle(1-p_{3}^{2})^{2}\rangle^{2}}=\frac{\mathcal{C}_{4}(Re_{\lambda})}{\mathcal{C}_{2}(Re_{\lambda})^{2}}$
(72)

to characterize the degree of non-Gaussianity in the PDF. From the inset of
figure 16(b), we can observe that the kurtosis monotonically increases with
$Re_{\lambda}$. Thus, the model PDF that we have developed in this section,
based on orientation dynamics in an ensemble of simple shear flows, captures
the non-Gaussian features observed in the experiments and appears promising in
exploring the suspension rheology of heavy anisotropic particles in
turbulence.

##  4 Conclusions

The model and the results presented in this paper are relevant in several
engineering and environmental scenarios involving sedimenting anisotropic
particles in turbulent flows. The dynamics of plankton in the oceans and
marine snow [14], fiber suspensions [41], pollen dispersion and particle or
cell aggregates in stirred tank reactors [15, 33] often involve competition
between the effects of gravitational settling and turbulence. The current
study focuses on the orientation dynamics of settling heavy particles in
homogeneous isotropic turbulence, accounting for fluid inertia due to
sedimentation and extending beyond single fiber models to more complex
ramified particles. The particles are assumed to be small enough, so angular
acceleration is negligible. This scenario is significant in atmospheric
research when studying the orientation distribution of ice crystals in cold
cirrus clouds ($T=-25^{\circ}$C to $-50^{\circ}$C). Heymsfield and Iaquinta
[25] presents an elaborate bullet rosette model (similar to ramified
particles) for different ice crystal shapes with maximum lengths ranging from
$L=100$ $\mu$m to $L=1$ mm. The corresponding terminal velocities for ice
crystals with aspect ratio $\kappa=20$ range from 2 cm/s ($L=100$ $\mu$m) to
80 cm/s ($L=1$ mm). This results in a particle Reynolds number of
$\textit{Re}_{\ell}=0.09$ up to $\textit{Re}_{\ell}=35$.

Estimates for turbulence intensities and energy dissipation rates for warm
cumulus clouds ($T=0^{\circ}$C to $10^{\circ}$C) can be found in Siebert et
al. [59] and Siebert et al. [60], with mean energy dissipation rates ranging
from $\langle\epsilon\rangle=\sim 10^{-3}$ m2/s3 to $\sim 10^{-2}$ m2/s3,
respectively, and with a fluid viscosity of about $\nu=1.4\times 10^{-5}$
m2/s). This yields Kolmogorov lengths of $\eta=1.1$ mm and $\eta=0.4$ mm,
which puts the ice crystals in the small particle limit, even with $L=2.5\eta$
Parsa and Voth [53]. The corresponding Kolmogorov velocities are $u_{\eta}=13$
mm/s and $u_{\eta}=32$ mm/s.

The ratio of ice crystal terminal velocities and Kolmogorov velocities enables
us to estimate the range of settling parameters $S_{F}$ for these atmospheric
conditions using Eq. 12. For the smallest ice crystal size, $S_{F}=0.4$
(Siebert et al. [59]) and $S_{F}=0.06$ (Siebert et al. [60]), and one could
expect no strong preferential alignment. However, the turbulence statistics
are taken from warm cumulus clouds, which are known to be more turbulent than
cirrus clouds. The settling parameter increases with particle size (terminal
velocity) to about $S_{F}\approx 600$ (Siebert et al. [59]) and $S_{F}\approx
100$ (Siebert et al. [60]) for the largest ice crystals. For all sizes except
the smallest ice crystals, $S_{F}\gg 1$ and we expect a strong alignment of
ice crystals with their preferential sedimentation orientation. For large ice
crystals, $S_{F}$ will be very large, and the orientation distributions are
most likely affected by particle asymmetries that prevent perfect alignment
rather than turbulence.

## Acknowledgments

This work was supported by the Army Research Office under grant
W911NF-15-1-0205. A.R. would like to acknowledge the support from Laboratory
for Atmospheric and Climate Sciences, Indian Institute of Technology Madras.

## References

  * Anand et al. [2020] Anand, P., Ray, S.S., Subramanian, G., 2020.  Orientation dynamics of sedimenting anisotropic particles in turbulence.  Physical Review Letters 125, 034501\.
  * Ardekani et al. [2017] Ardekani, M.N., Sardina, G., Brandt, L., Karp-Boss, L., Bearon, R., Variano, E., 2017\.  Sedimentation of inertia-less prolate spheroids in homogenous isotropic turbulence with application to non-motile phytoplankton.  Journal of Fluid Mechanics 831, 655–674.
  * Batchelor [1970] Batchelor, G.K., 1970.  Slender-body theory for particles of arbitrary cross-section in stokes flow.  Journal of Fluid Mechanics 44, 419–440.
  * Bosse et al. [2006] Bosse, T., Kleiser, L., Meiburg, E., 2006.  Small particles in homogeneous turbulence: Settling velocity enhancement by two-way coupling.  Physics of Fluids 18, 027102\.
  * Bragg et al. [1974] Bragg, G.M., van Zuiden, L., Hermance, C.E., 1974.  The free fall of cylinders at intermediate reynold’s numbers.  Atmospheric Environment (1967) 8, 755–764.  URL: http://www.sciencedirect.com/science/article/pii/0004698174901656, doi:http://dx.doi.org/10.1016/0004-6981(74)90165-6.
  * Brenner [1963] Brenner, H., 1963.  The stokes resistance of an arbitrary particle.  Chemical Engineering Science 18, 1–25.
  * Brunk et al. [1998] Brunk, B.K., Koch, D.L., Lion, L.W., 1998.  Turbulent coagulation of colloidal particles.  Journal of Fluid Mechanics 364, 81–113.
  * Challabotla et al. [2015] Challabotla, N.R., Zhao, L., Andersson, H.I., 2015.  Orientation and rotation of inertial disk particles in wall turbulence.  Journal of Fluid Mechanics 766.  doi:10.1017/jfm.2015.38.
  * Cho et al. [1981] Cho, H.R., Iribarne, J.V., Richards, W.G., 1981.  On the orientation of ice crystals in a cumulonimbus cloud.  Journal of the Atmospheric Sciences 38, 1111–1114.  URL: http://dx.doi.org/10.1175/1520-0469(1981)038<1111:OTOOIC>2.0.CO;2, doi:10.1175/1520-0469(1981)038<1111:OTOOIC>2.0.CO;2.
  * Cox [1965] Cox, R.G., 1965.  The steady motion of a particle of arbitrary shape at small reynolds numbers.  Journal of Fluid Mechanics 23, 625–643.  doi:10.1017/S0022112065001593.
  * Dabade et al. [2015] Dabade, V., Marath, N.K., Subramanian, G., 2015.  Effects of inertia and viscoelasticity on sedimenting anisotropic particles.  Journal of Fluid Mechanics 778, 133–188.
  * Devenish et al. [2012] Devenish, B.J., Bartello, P., Brenguier, J.L., Collins, L.R., Grabowski, W.W., IJzermans, R.H.A., Malinowski, S.P., Reeks, M.W., Vassilicos, J.C., Wang, L.P., Warhaft, Z., 2012.  Droplet growth in warm turbulent clouds.  Quarterly Journal of the Royal Meteorological Society 138, 1401–1429.  URL: http://dx.doi.org/10.1002/qj.1897, doi:10.1002/qj.1897.
  * Dhanasekaran et al. [2021] Dhanasekaran, J., Roy, A., Koch, D.L., 2021.  Collision rate of bidisperse, hydrodynamically interacting spheres settling in a turbulent flow.  Journal of Fluid Mechanics 912, A5.
  * Dunlop et al. [1994] Dunlop, E.H., Namdev, P.K., Rosenberg, M.Z., 1994.  Effect of fluid shear forces on plant cell suspensions.  Chemical Engineering Science 49, 2263–2276.
  * Ehrl et al. [2008] Ehrl, L., Soos, M., Morbidelli, M., 2008.  Dependence of aggregate strength, structure, and light scattering properties on primary particle size under turbulent conditions in stirred tank.  Langmuir 24, 3070–3081.
  * Ern et al. [2012] Ern, P., Risso, F., Fabre, D., Magnaudet, J., 2012.  Wake-induced oscillatory paths of bodies freely rising or falling in fluids.  Annual Review of Fluid Mechanics 44, 97–121.
  * Fang et al. [2016] Fang, L., Zhang, Y., Fang, J., Zhu, Y., 2016.  Relation of the fourth-order statistical invariants of velocity gradient tensor in isotropic turbulence.  Physical Review E 94, 023114\.
  * Field et al. [1997] Field, S.B., Klaus, M., Moore, M.G., Nori, F., 1997\.  Chaotic dynamics of falling disks.  Nature 388, 252–254.
  * Girimaji and Pope [1990] Girimaji, S., Pope, S., 1990\.  A diffusion model for velocity gradients in turbulence.  Physics of Fluids A: Fluid Dynamics 2, 242–256.
  * Gustavsson et al. [2017] Gustavsson, K., Jucha, J., Naso, A., Lévêque, E., Pumir, A., Mehlig, B., 2017\.  Statistical model for the orientation of nonspherical particles settling in turbulence.  Physical review letters 119, 254501\.
  * Gustavsson et al. [2019] Gustavsson, K., Sheikh, M., Lopez, D., Naso, A., Pumir, A., Mehlig, B., 2019\.  Effect of fluid inertia on the orientation of a small prolate spheroid settling in turbulence.  New Journal of Physics 21, 083008\.
  * Hashino et al. [2014] Hashino, T., Chiruta, M., Polzin, D., Kubicek, A., Wang, P.K., 2014.  Numerical simulation of the flow fields around falling ice crystals with inclined orientation and the hydrodynamic torque.  Atmospheric Research 150, 79–96.  URL: http://www.sciencedirect.com/science/article/pii/S0169809514002658, doi:https://doi.org/10.1016/j.atmosres.2014.07.003.
  * Hejazi et al. [2017] Hejazi, B., Mehlig, B., Voth, G.A., 2017.  Emergent scar lines in chaotic advection of passive directors.  arXiv preprint arXiv:1706.07398v1 .
  * Henthorn et al. [2005] Henthorn, K.H., Park, K., Curtis, J.S., 2005.  Measurement and prediction of pressure drop in pneumatic conveying: Effect of particle characteristics, mass loading, and reynolds number.  Industrial & engineering chemistry research 44, 5090–5098.
  * Heymsfield and Iaquinta [2000] Heymsfield, A.J., Iaquinta, J., 2000\.  Cirrus crystal terminal velocities.  Journal of the atmospheric sciences 57, 916–938.
  * Hu [2007] Hu, Y., 2007.  Depolarization ratio–effective lidar ratio relation: Theoretical basis for space lidar cloud phase discrimination.  Geophysical Research Letters 34.
  * Ireland et al. [2016] Ireland, P.J., Bragg, A.D., Collins, L.R., 2016.  The effect of reynolds number on inertial particle dynamics in isotropic turbulence. part 2. simulations with gravitational effects.  Journal of Fluid Mechanics 796, 659–711.
  * Jayaweera and Mason [1965] Jayaweera, K.O.L.F., Mason, B.J., 1965\.  The behaviour of freely falling cylinders and cones in a viscous fluid.  Journal of Fluid Mechanics 22, 709–720.  doi:10.1017/S002211206500109X.
  * Jeffery [1922] Jeffery, G.B., 1922.  The motion of ellipsoidal particles immersed in a viscous fluid.  Proceedings of the Royal Society of London. Series A, Containing papers of a mathematical and physical character 102, 161–179.
  * Kajikawa [1992] Kajikawa, M., 1992.  Observations of the falling motion of plate-like snow crystals part i: The free-fall patterns and velocity.  Journal of the Meteorological Society of Japan. Ser. II 70, 1–9.
  * Khayat and Cox [1989] Khayat, R.E., Cox, R.G., 1989\.  Inertia effects on the motion of long slender bodies.  Journal of Fluid Mechanics 209, 435–462.  doi:10.1017/S0022112089003174.
  * Kim and Karrila [2013] Kim, S., Karrila, S.J., 2013\.  Microhydrodynamics: principles and selected applications.  Courier Corporation.
  * Kiørboe [2001] Kiørboe, T., 2001.  Formation and fate of marine snow: small-scale processes with large-scale implications.  Scientia marina 65, 57–71.
  * Klett [1995] Klett, J.D., 1995.  Orientation model for particles in turbulence.  Journal of the Atmospheric Sciences 52, 2276–2285.
  * Koch and Pope [2002] Koch, D.L., Pope, S.B., 2002\.  Coagulation-induced particle-concentration fluctuations in homogeneous, isotropic turbulence.  Physics of Fluids 14, 2447–2455.
  * Kramel [2017] Kramel, S., 2017.  Non-spherical particle dynamics in turbulence.  Ph.D. thesis. Wesleyan University. https://digitalcollections.wesleyan.edu/object/ir-2304.
  * Kramel et al. [2016] Kramel, S., Voth, G.A., Tympel, S., Toschi, F., 2016\.  Preferential rotation of chiral dipoles in isotropic turbulence.  Physical Review Letters 117, 154501\.  URL: https://link.aps.org/doi/10.1103/PhysRevLett.117.154501, doi:10.1103/PhysRevLett.117.154501.
  * Krushkal and Gallily [1988] Krushkal, E., Gallily, I., 1988\.  On the orientation distribution function of nonspherical aerosol particles in a general shear flow–ii. the turbulent case.  Journal of Aerosol Science 19, 197–211.
  * Leal [1980] Leal, L.G., 1980.  Particle motions in a viscous fluid.  Annual Review of Fluid Mechanics 12, 435–476.  URL: http://dx.doi.org/10.1146/annurev.fl.12.010180.002251, doi:10.1146/annurev.fl.12.010180.002251.
  * Lopez and Guazzelli [2017] Lopez, D., Guazzelli, É., 2017\.  Inertial effects on fibers settling in a vortical flow.  Physical Review Fluids 2, 024306\.
  * Lundell et al. [2011] Lundell, F., Söderberg, L.D., Alfredsson, P.H., 2011.  Fluid mechanics of papermaking.  Annual Review of Fluid Mechanics 43, 195–217.
  * Marchioli et al. [2010] Marchioli, C., Fantoni, M., Soldati, A., 2010.  Orientation, distribution, and deposition of elongated, inertial fibers in turbulent channel flow.  Physics of Fluids 22, 033301\.  URL: http://dx.doi.org/10.1063/1.3328874, doi:10.1063/1.3328874.
  * Marcus et al. [2014] Marcus, G.G., Parsa, S., Kramel, S., Ni, R., Voth, G.A., 2014.  Measurements of the solid-body rotation of anisotropic particles in 3d turbulence.  New Journal of Physics 16, 102001\.
  * Meneveau [2011] Meneveau, C., 2011.  Lagrangian dynamics and models of the velocity gradient tensor in turbulent flows.  Annual Review of Fluid Mechanics 43, 219–245.
  * Menon et al. [2017] Menon, U., Roy, A., Kramel, S., Voth, G., Koch, D., Collaboration, V.L., et al., 2017\.  Theoretical predictions of the orientation distribution of high-aspect-ratio, inertial particles settling in isotropic turbulence, in: APS Division of Fluid Dynamics Meeting Abstracts, pp. Q36–011.
  * Menon [2019] Menon, U.K., 2019.  Theory and Simulation for the Orientation of High Aspect Ratio Particles Settling in Homogeneous Isotropic Turbulence.  Master’s thesis. Cornell University. https://doi.org/10.7298/8kfe-6s23.
  * Mortensen et al. [2008] Mortensen, P.H., Andersson, H.I., Gillissen, J.J.J., Boersma, B.J., 2008\.  On the orientation of ellipsoidal particles in a turbulent shear flow.  International Journal of Multiphase Flow 34, 678–683.  URL: http://www.sciencedirect.com/science/article/pii/S0301932208000098, doi:http://dx.doi.org/10.1016/j.ijmultiphaseflow.2007.12.007.
  * Newsom and Bruce [1994] Newsom, R.K., Bruce, C.W., 1994\.  The dynamics of fibrous aerosols in a quiescent atmosphere.  Physics of Fluids 6, 521–530.  URL: http://dx.doi.org/10.1063/1.868347, doi:10.1063/1.868347.
  * Newsom and Bruce [1998] Newsom, R.K., Bruce, C.W., 1998\.  Orientational properties of fibrous aerosols in atmospheric turbulence.  Journal of Aerosol Science 29, 773–797.
  * Noel et al. [2002] Noel, V., Roy, G., Bissonnette, L., Chepfer, H., Flamant, P., 2002.  Analysis of lidar measurements of ice clouds at multiple incidence angles.  Geophysical Research Letters 29.
  * Noel and Sassen [2005] Noel, V., Sassen, K., 2005.  Study of planar ice crystal orientations in ice clouds from scanning polarization lidar observations.  Journal of Applied Meteorology 44, 653–664.
  * Parsa et al. [2012] Parsa, S., Calzavarini, E., Toschi, F., Voth, G.A., 2012\.  Rotation rate of rods in turbulent fluid flow.  Physical Review Letters 109, 134501\.
  * Parsa and Voth [2014] Parsa, S., Voth, G.A., 2014\.  Inertial range scaling in rotations of long rods in turbulence.  Physical Review Letters 112, 024501\.
  * Pumir and Wilkinson [2011] Pumir, A., Wilkinson, M., 2011\.  Orientation statistics of small particles in turbulence.  New Journal of Physics 13, 093030\.
  * Roy et al. [2019] Roy, A., Hamati, R.J., Tierney, L., Koch, D.L., Voth, G.A., 2019.  Inertial torques and a symmetry breaking orientational transition in the sedimentation of slender fibres.  Journal of Fluid Mechanics 875, 576–596.
  * Sabban et al. [2017] Sabban, L., Cohen, A., van Hout, R., 2017.  Temporally resolved measurements of heavy, rigid fibre translation and rotation in nearly homogeneous isotropic turbulence.  Journal of Fluid Mechanics 814, 42–68.
  * Shin and Koch [2005] Shin, M., Koch, D.L., 2005.  Rotational and translational dispersion of fibres in isotropic turbulent flows.  Journal of Fluid Mechanics 540, 143–173.
  * Shin et al. [2006] Shin, M., Koch, D.L., Subramanian, G., 2006.  A pseudospectral method to evaluate the fluid velocity produced by an array of translating slender fibers.  Physics of Fluids 18, 063301\.  URL: http://aip.scitation.org/doi/abs/10.1063/1.2205200, doi:10.1063/1.2205200.
  * Siebert et al. [2006] Siebert, H., Lehmann, K., Wendisch, M., Shaw, R., 2006\.  Small-scale turbulence in clouds, in: 12th Conference on Cloud Physics, pp. 10–14.
  * Siebert et al. [2015] Siebert, H., Shaw, R., Ditas, J., Schmeissner, T., Malinowski, S., Bodenschatz, E., Xu, H., 2015.  High-resolution measurement of cloud microphysics and turbulence at a mountaintop station.  Atmospheric Measurement Techniques 8, 3219–3228.
  * Siewert et al. [2014] Siewert, C., Kunnen, R., Meinke, M., Schröder, W., 2014\.  Orientation statistics and settling velocity of ellipsoids in decaying turbulence.  Atmospheric research 142, 45–56.
  * Siggia [1981] Siggia, E.D., 1981.  Invariants for the one-point vorticity and strain rate correlation functions.  The Physics of Fluids 24, 1934–1936.
  * Vakil and Green [2009] Vakil, A., Green, S.I., 2009\.  Drag and lift coefficients of inclined finite circular cylinders at moderate reynolds numbers.  Computers & Fluids 38, 1771–1781.
  * Variano and Cowen [2008] Variano, E.A., Cowen, E.A., 2008\.  A random-jet-stirred turbulence tank.  Journal of Fluid Mechanics 604, 1–32.
  * Voth and Soldati [2017] Voth, G.A., Soldati, A., 2017\.  Anisotropic particles in turbulence.  Annual Review of Fluid Mechanics 49, null.  URL: http://dx.doi.org/10.1146/annurev-fluid-010816-060135, doi:10.1146/annurev-fluid-010816-060135.
  * Wang and Maxey [1993] Wang, L.P., Maxey, M.R., 1993\.  Settling velocity and concentration distribution of heavy particles in homogeneous isotropic turbulence.  Journal of Fluid Mechanics 256, 27–68.
  * Westbrook [2011] Westbrook, C.D., 2011.  Origin of the parry arc.  Quarterly Journal of the Royal Meteorological Society 137, 538–543.
  * Westbrook et al. [2010] Westbrook, C.D., Illingworth, A.J., O’Connor, E.J., Hogan, R.J., 2010\.  Doppler lidar measurements of oriented planar ice crystals falling from supercooled and glaciated layer clouds.  Quarterly Journal of the Royal Meteorological Society 136, 260–276.
  * Wilkinson et al. [2009] Wilkinson, M., Bezuglyy, V., Mehlig, B., 2009.  Fingerprints of random flows.  Physics of Fluids 21, 043304\.
  * Willmarth et al. [1964] Willmarth, W.W., Hawk, N.E., Harvey, R.L., 1964.  Steady and unsteady motions and wakes of freely falling disks.  Physics of Fluids 7, 197–208.
  * Yeung and Pope [1989] Yeung, P., Pope, S., 1989.  Lagrangian statistics from direct numerical simulations of isotropic turbulence.  Journal of Fluid Mechanics 207, 531–586.
  * Zastawny et al. [2012] Zastawny, M., Mallouppas, G., Zhao, F., van Wachem, B., 2012\.  Derivation of drag and lift force and torque coefficients for non-spherical particles in flows.  International Journal of Multiphase Flow 39, 227–239.  URL: http://www.sciencedirect.com/science/article/pii/S0301932211002047, doi:https://doi.org/10.1016/j.ijmultiphaseflow.2011.09.004.
  * Zhao and Andersson [2016] Zhao, L., Andersson, H.I., 2016\.  Why spheroids orient preferentially in near-wall turbulence.  Journal of Fluid Mechanics 807, 221–234.
  * Zhao et al. [2010] Zhao, L.H., Andersson, H.I., Gillissen, J.J.J., 2010.  Turbulence modulation and drag reduction by spherical particles.  Physics of Fluids 22, 081702\.  URL: http://dx.doi.org/10.1063/1.3478308, doi:10.1063/1.3478308.
  * Zikmunda and Vali [1972] Zikmunda, J., Vali, G., 1972\.  Fall patterns and fall velocities of rimed ice crystals.  Journal of the Atmospheric Sciences 29, 1334–1347.