# Near-wall turbulence modulation by small inertial particles

Pedro Costa1<EMAIL_ADDRESS>Luca Brandt2  Francesco Picano3
1Faculty of Industrial Engineering, Mechanical Engineering and Computer
Science, University of Iceland, Hjardarhagi 2-6, 107 Reykjavik, Iceland
2Linné FLOW Centre and SeRC (Swedish e-Science Research Centre), Department of
Engineering Mechanics, KTH, SE-100 44 Stockholm, Sweden
3Department of Industrial Engineering, University of Padova, Via Venezia, 1,
35131, Padova, Italy

###### Abstract

We use interface-resolved simulations to study near-wall turbulence modulation
by small inertial particles, much denser than the fluid, in dilute/semi-dilute
conditions. We considered three bulk solid mass fractions, $\Psi=0.34\%$,
$3.37\%$ and $33.7\%$, with only the latter two showing turbulence modulation.
The increase of the drag is strong at $\Psi=3.37\%$, but mild in the densest
case. Two distinct regimes of turbulence modulation emerge: for smaller mass
fractions, the turbulence statistics are weakly affected and the near-wall
particle accumulation increases the drag so the flow appears as a single-phase
flow at slightly higher Reynolds number. Conversely, at higher mass fractions,
the particles modulate the turbulent dynamics over the entire flow, and the
interphase coupling becomes more complex. In this case, fluid Reynolds
stresses are attenuated, but the inertial particle dynamics near the wall
increases the drag via correlated velocity fluctuations, leading to an overall
drag increase. Hence, we conclude that, although particles at high mass
fractions reduce the fluid turbulent drag, the solid phase inertial dynamics
still increases the overall drag. However, inspection of the streamwise
momentum budget in the two-way coupling limit of vanishing volume fraction,
but finite mass fraction, indicates that this trend could reverse at even
higher particle load.

##  1 Introduction

The turbulent channel flow laden with small, gravity-free inertial spheres is
a paradigmatic multiphase flow, yet not fully understood and correctly
modelled. This system has been widely studied using direct numerical
simulations (DNS) of the Navier-Stokes equations using the point-particle
approximation (Balachandar & Eaton, 2010), which assumes that interphase
coupling is localized at a single point and the particle dynamics driven by an
undisturbed flow field, evaluated at the particle position. This approximation
allows us to study the flow dynamics beyond the particle scale, without
solving the actual flow around each particle. Despite the many studies, the
validity of these particle–fluid coupling models has not been fully assessed
for canonical turbulent wall flows, in particular when the underlying
turbulence is altered by the dispersed phase.

When the particle feedback on the flow becomes important – the two-way
coupling regime – care should be taken in the estimation of the undisturbed
velocity. The challenge is conciliating the estimation of an undisturbed
velocity sampled by the particle with the need of forcing the local velocity
field. Indeed, approaches for accurate and realistic two-way coupling methods
have been pursued since the work of Crowe et al. (1977), and are still object
of active research (see, e.g., Gualtieri et al., 2015; Horwitz & Mani, 2016;
Ireland & Desjardins, 2017; Battista et al., 2019; Pakseresht et al., 2020).
Moreover, results in the literature for integral observables show
qualitatively different trends. For instance, while Vreman (2007); Zhao et al.
(2010) have observed a drag-reducing behaviour in two-way coupled particle-
laden turbulent channel flow, other studies have measured a drag increase
(Battista et al., 2019). These opposite trends suggest the presence of
different regimes of momentum transfer in wall-bounded particle-laden
turbulence. Indeed, Capecelatro et al. (2018) used Euler–Lagrange volume-
filtered DNS to investigate in detail the turbulent kinetic energy (TKE)
budget in vertical particle-laden channel flow at different volume fractions.
Different regimes were observed with increasing mass fraction – first, a
turbulent regime at low mass fractions; second, increasing the mass fraction,
the turbulence activity decreases until the flow laminarizes, due to the
growing importance of a “drag dissipation-and-exchange term” associated with
fluid–particle correlated velocity fluctuations; finally, increasing the mass
fraction beyond this limit triggers a second mechanism for TKE increase, due
to the average velocity difference between phases, which grows important and
re-energizes fluid turbulence. Despite the exciting developments exploiting
particle-modelled two-way coupling interactions, there is a clear need for
high-fidelity data to support the body of work using these modelled
simulations, and to possibly reconcile seemingly conflicting observations of
drag increase and drag reduction in two-way coupling point-particle numerical
simulations.

Indeed, the tremendous developments of approaches for DNS with modelled
particle–fluid coupling have not been accompanied by the high-fidelity data
which are essential for validating the underlying assumptions. These data can
be obtained either from interface-resolved simulations of the Navier-Stokes
equations, or well-controlled experiments matching the computational set-up.
Nowadays, the first direct comparisons between point-particle DNS data and
experiments or interface-resolved DNS are possible. Examples are the work by
Wang et al. (2019), which compares one-to-one two-way coupled point-particle
DNS with experimental measurements for vertical turbulent channel flow at
different volume loads, and the interface-resolved (also denoted particle-
resolved) DNS of spherical particles near the point-particle limit in
homogeneous isotropic turbulence (Schneiders et al., 2017; Mehrabadi et al.,
2018; Fröhlich et al., 2018) and in turbulent channel flow (Horne & Mahesh,
2019; Costa et al., 2020).

The present work also exploits interface-resolved simulations and focuses on
near-wall turbulence modulation by small inertial particles. We adopt the
configuration in Costa et al. (2020) and add results at higher volume
fraction. Three bulk mass fractions $\Psi=0.34\%$, $3.37\%$ and $33.7\%$ at
fixed particle to fluid density ratio ($\Pi_{\rho}=100$) are hence considered,
where the latter two show non-negligible two-way coupling effects. Our results
reveal two distinct mechanisms for turbulence modulation. At lower volume
fractions, the flow is only dense very close to the wall, and the higher flow
inertia in this region, with particles travelling much faster than the fluid,
results in a drag-increased flow resembling single-phase turbulence at
slightly higher Reynolds number. At higher volume fractions, the dispersed
phase is dynamically important over the entire channel. Here, the drag
increasing effect is amplified by correlated particle velocity fluctuations,
but counterbalanced by a substantial fluid turbulence drag reduction. This
results in a milder increase in drag for the denser case. These observations
may be explained in light of the streamwise momentum balance for vanishing
volume fraction, but with non-negligible mass fraction.

##  2 Methods and Computational Set-up

Since we use here the tools and set-up in Costa et al. (2020), with one
additional case at the largest volume fraction, we briefly summarize the
numerical method and refer to the previous work for more details. We solve the
continuity and Navier-Stokes equations for an incompressible Newtonian fluid,
together with the Newton-Euler equations driving the motion of the solid
spherical particles. These two sets of equations are coupled directly using
the immersed-boundary method developed by Breugem (2012) (see also Uhlmann,
2005), built on a standard second-order finite-difference method on a three-
dimensional, staggered Cartesian grid, using a fast-Fourier-transform-based
pressure-projection method (Kim & Moin, 1985; Costa, 2018). Short-range
particle–particle/particle–wall interactions (lubrication and solid–solid
contact) are modelled using the method of Costa et al. (2015), as in Costa et
al. (2020). More specifically, as regards solid–solid contact, the particles
are frictionless, and with a normal solid–solid coefficient of restitution of
$0.97$.

Turbulent channel flow is simulated in a domain periodic in the streamwise
($x$) and spanwise ($z$) directions, with no-slip and no-penetration boundary
conditions imposed at the walls ($y=h\mp h$), where $h$ is the channel half-
height. The flow is driven by a uniform pressure gradient that ensures a
constant bulk velocity. The bulk Reynolds number is equal to
$\Rey_{b}=U_{b}(2h)/\nu=5\,600$, which corresponds to an unladen friction
Reynolds number $\Rey_{\tau}^{sph}=u_{\tau}h/\nu\approx 180$; where $U_{b}$ is
the bulk flow velocity and $u_{\tau}$ the wall friction velocity. The particle
properties are chosen to yield a particle Reynolds number
$\Rey_{p}=Du_{\tau}/\nu=D^{+}=3$, and Stokes number
$St_{p}=\Pi_{\rho}\Rey_{p}^{2}/18=50$, where $\nu$ is the fluid kinematic
viscosity, $D$ is the particle diameter and $\Pi_{\rho}$ the particle-to-fluid
mass density ratio; since the particle size is restricted by resolution
requirements, $\Pi_{\rho}$ is used to achieve the target particle Stokes
number, which was demonstrated to feature highly inhomogeneous particle
distributions in wall turbulence (see, e.g., Sardina et al., 2012). Three
values of solid volume fraction are considered, increasing by factors of $10$:
$\Phi\simeq 3\cdot 10^{-5}$, denoted very dilute (VD), $\Phi\simeq 3\cdot
10^{-4}$, dilute (D), and $\Phi\simeq 3\cdot 10^{-3}$, semi-dilute (SD), with
the data pertaining to the first two cases also used in the recent study by
Costa et al. (2020). Table 1 shows all relevant physical and computational
parameters.

Case | $\Phi\,\,\,\,(N_{p})$ | $\Psi$ | Notes
---|---|---|---
VD | $0.003\%\,\,(500)\hphantom{0}\hphantom{0}\hphantom{0}$ | $0.337\%$ | interface-resolved (very dilute)
D0 | $0.034\%\,\,(5\,000)\hphantom{0}\hphantom{0}$ | $3.367\%$ | interface-resolved (dilute)
SD | $0.337\%\,\,(50\,000)\hphantom{0}$ | $33.67\%$ | interface-resolved (semi dilute)
Table 1: Computational parameters. $\Phi$/$\Psi$ denote the bulk solid
volume/mass fraction, and $N_{p}$ the total number of particles. Common to all
cases: bulk Reynolds number $\Rey_{b}=5\,600$ (i.e. friction Reynolds number
in the single-phase limit $\Rey_{\tau}^{sph}\approx 180$); particle size ratio
$D/(2h)=1/120$; particle-to-fluid mass density ratio $\Pi_{\rho}=100$. These
correspond to a particle Reynolds number based on the unladen reference values
of $D^{+}=3$ and Stokes number $St=50$. The fluid domain is discretized on a
regular Cartesian grid with
$(L_{x}/N_{x})\times(L_{y}/N_{y})\times(L_{z}/N_{z})=(6h/4320)\times(2h/1440)\times(3h/2160)$,
while the particles are resolved with $D/\Delta x=12$ grid points over the
particle diameter (in total $420$ Lagrangian grid points) (same as Costa et
al., 2020).

Finally, unless otherwise stated, the mesoscale-averaged profiles reported in
this manuscript correspond to intrinsic averages (i.e., averaged only over the
corresponding phase) in time and along the two homogeneous directions for each
phase (Costa et al., 2020), obtained with wall-parallel bins with thickness
equal to the grid spacing.

##  3 Results

A three-dimensional visualization of the flow pertaining to the different
cases is reported in figure 1, where iso-surfaces of positive second invariant
of the velocity gradient tensor $Q$ (coloured by the local wall-normal
velocity) and the particles are displayed. High-vorticity spots, footprints of
the presence of the particles, are noticeable in all particle-laden cases,
especially closer to the wall, where, as we will see, the mean slip velocity
between the phases is highest. As expected, no significant qualitative
differences can be seen between the single-phase and the very dilute case VD,
apart from the very small number of dispersed particles in the latter case.
Qualitative differences between these cases and case D are also small;
however, an increased number of high-vorticity spots due to the particles is
obvious over the entire domain. Finally, a strong modulation of the flow
dynamics by the particles is seen when inspecting case SD, where the
disruption of the flow coherent structures is evident. One of the main
messages of the present work is that, although cases D and SD have macroscopic
flow variations and, therefore, are formally in a two-way coupling regime, the
mechanism of turbulence modulation is definitively different.

single-phasecase VD
case  Dcase SD

Figure 1: Three-dimensional flow visualizations of the different cases under
consideration. We display surfaces of constant second invariant of the
velocity gradient tensor $Q=20(U_{b}/h)^{2}$, coloured by the local wall-
normal velocity (white – high and blue – low). The particles are shown to
scale, in orange colour. See table 1 for a description of the different cases.

Figure 2(a) presents the friction Reynolds number as a function of the bulk
volume fraction, while the change in drag relative to the unladen flow is
quantified in panel (b). The friction Reynolds number increases significantly
($5\%$) from case VD to case D, despite the relatively low mass fraction
$\Psi=\mathcal{O}(10^{-2})$. As discussed in Costa et al. (2020), this
increase is attributed to the higher inertia induced by the large local mass
fraction near the wall; the particles are driven by turbophoresis towards the
wall, where they experience a large apparent slip velocity. Further increasing
the mass fraction by an order of magnitude (case SD) results in a mild
increase in drag, only about $3\%$ with respect to case D. This milder
increase indicates that a competing drag-reducing mechanism comes into play at
higher mass loading. The blue symbols in figure 2 display quantities computed
using a wall friction velocity determined from the centreline slope of the
Reynolds stress, so as to quantify the decrease in turbulent momentum
transfer, as discussed in detail later.

(a)(b)

Figure 2: (a) Friction Reynolds number $\Rey_{\tau}$ versus the bulk volume
fraction $\Psi$; $\Rey_{\tau}^{t}$ is based on a velocity scale computed from
the slope of the outer-scaled Reynolds shear stresses evaluated at the channel
centreline $u_{\tau}^{t}=\sqrt{\partial_{y/h}\left\langle
u^{\prime}v^{\prime}\right\rangle|_{y/h=1}}$, as a rough measure of the
contribution of fluid velocity correlations to the mean wall shear (see Picano
et al., 2015). (b) Drag modulation
$\mathrm{DC}=\Rey_{\tau}^{2}/\Rey_{\tau}^{2,sph}-1$, computed with
$\Rey_{\tau}$ and $\Rey_{\tau}^{t}$ ($\mathrm{DC}^{t}$).

Figure 3(a) shows the mass fraction normalized by the corresponding bulk value
versus the wall-normal distance in particle diameters. All the cases show a
peak at the wall, corresponding to a particle layer. Note also that these
peaks occur slightly beyond a wall distance of one particle radius, due to the
bouncing dynamics caused by particle–wall collisions. Clearly, the fraction of
particles near the wall decreases with increasing mass loading. This decrease
in wall accumulation is attributed to the increased shear rate near the
boundary, which enhances lift forces (Costa et al., 2020). However, case SD
shows a disproportionally strong decrease of the relative particle
concentration, considering the relatively low increase in wall shear with
respect to case D. This is partly caused by two-way coupling effects, which
dampen the intensity of the wall-normal fluid velocity fluctuations
responsible – to first approximation – for the turbophoretic drift (Marchioli
& Soldati, 2002). In addition, particle–particle interactions may become
significant. The inset of figure 3(a) shows the mass fraction profile as a
function of the outer-scaled wall-normal distance. While the mass fraction
profiles become uniform far from the wall for cases VD and D, the distribution
shows a mild monotonic increase with the wall distance in case SD. These
observations indicate that, for this densest case, non-negligible
particle–particle interactions may be driving particles towards regions of low
shear. The mechanism underlying this drift is most likely similar to that
reported in Fornari et al. (2016) – particle inertia and local high shear
promote inter-particle interactions, causing a net migration towards low shear
regions. We attribute the non-zero slope at the channel centreline of the
profile pertaining to case SD to a lower statistical convergence in the
channel core in combination with the amplification induced by the logarithmic
scale, as explained in the figure caption.

(a)(b)

Figure 3: (a) Local solid mass fraction $\psi$, normalized by the bulk value
$\Psi$ as a function of the wall-normal distance in particle diameters (note
that the profile is clamped beyond $y/D=6$, to highlight the near-wall peaks);
the inset shows the nominal (i.e., not normalized) local mass fraction $\psi$
as a function of the outer-scaled wall-normal distance. The non-zero slope at
the centreline, for case SD, is attributed to insufficient statistical
convergence in this region where the flow is more dilute and the volume
fraction fluctuations higher, further highlighted by the logarithmic scaling.
(b) Outer-scaled spanwise autocorrelation of streamwise velocity at
$y/h=0.067$ and $y/h=0.167$, ($y/\delta_{\nu}^{sph}\approx 12$ and $30$);
‘sph’ corresponds to the unladen flow.

The inset of figure 3(a) also illustrates where two-way coupling effects are
expected to be important. These can be envisaged near the wall for the less
dilute cases, as the solid mass fraction increases to $20\%$ for case D, and
to about $80\%$ for case SD. Away from the wall, instead, the solid mass
fraction retains high values only for case SD, $25\%-40\%$, and therefore
turbulence modulation is also expected in the bulk, which we will relate to
the milder overall drag increase from case D to case SD. In the latter case,
we observe a significant turbulence modulation by the solid particles: the
spacing between low- and high-speed streaks is larger, the spanwise velocity
variation smoother, and the maximum streak amplitude occurs at a much larger
wall-normal distance, which is reflected in the spanwise autocorrelation of
the streamwise velocity in wall-parallel planes shown in figure 3(b). This
trend can be clearly seen in figure 4, showing the contours of the near-wall
streamwise fluid velocity.

single-phasecase VDcase Dcase SD$u$

Figure 4: Contours of outer-scaled streamwise fluid velocity $u$ at a plane
$y/h=0.07$, together with the particles located at $y/h<0.15$. Perspective of
an observer aligned with the negative wall-normal direction, with flow from
left to right (the contour is shown with slight transparency).

Case SD also presents non-negligible four-way coupling effects. In particular,
we estimated the particle collision frequency near the wall following Sundaram
& Collins (1997):

$N_{c}=\frac{1}{2}\pi R^{2}n^{2}g(R)\left\langle\Delta
U^{n-}\right\rangle\mathrm{,}$ (1)

where $R$ is the particle radius, $n$ the local number density, $g(R)$ the
particle radial distribution function at the imminence of contact and
$\left\langle\Delta U^{n-}\right\rangle$ the mean particle relative velocity
projected along the line-of-centres, conditioned to negative values (i.e.
promoting a collision) at the imminence of contact (i.e. interparticle
distance equal to one particle diameter). We found the collision frequency for
case SD to be $5.5$ near the wall and to reduce to $0.53$ in the bulk, a
quantity defined per unit volume $L_{x}L_{z}D$ and bulk eddy turnover time
$(O(h/u_{\tau}))$, more than one order of magnitude larger than in case D
($0.29$ near the wall and $0.002$ in the bulk). We should note that, owing to
the computational cost and the low volume fractions considered, we were unable
to collect enough data to directly measure particle collision statistics.
Also, a direct assessment of collision frequency can yield deceiving results
due to the possibility of sustained contacts; see Kuerten & Vreman (2016).

The inner-scaled mean velocity profiles, see figure 5(a), resemble those of
the single-phase flow, with a downward shift that indicates drag increase (the
difference is less apparent in the corresponding outer-scaled profiles shown
in panel (b)). The inset of the figure shows the particle-to-fluid apparent
slip velocity, defined as the difference between the mean velocity profiles of
the fluid and solid phase. The negative minimum in the near-wall region is due
to the higher particle velocity where the fluid velocity is vanishing. The
local maximum of positive slip velocity occurs in the buffer region, which is
an indirect indication of the well-known tendency of near-wall inertial
particles to oversample regions of low streamwise fluid velocity observed in
numerous studies (e.g., Rouson & Eaton, 2001; Marchioli & Soldati, 2002, among
others). The higher the mass loading, the weaker this slip velocity is. Case
SD, in particular, shows a deeper minimum of the slip velocity, and virtually
no slip between the two phases beyond $y>10\delta_{\nu}$. This deeper minimum
is reflected in more pronounced wakes for particles moving faster than the
fluid, see figure 4. Interestingly, the fluid-to-particle apparent slip
velocity does not attain positive values for case SD, meaning that particles
are less prone to sampling low-speed regions. Since the low (and high) speed
regions are located further away from the wall, particles departing from the
wall may not reside for long enough in these regions. Finally, we should note
that wider streak spacing with larger wall-normal extent is often a feature of
turbulent drag reduction (see e.g. Tiederman et al., 1985).

(a)(b)

Figure 5: Inner- (a) and outer-scaled (b) mean streamwise velocity profiles.
Lines – fluid velocity; symbols – particle velocity. The insets show the
apparent fluid-to-particle slip velocity, defined as the difference between
the intrinsic mean velocity of each phase.

(a)(b)
(c)(d)

Figure 6: Second-order moments of mean fluid (lines) and particle (symbols)
velocity. (a) Inner-scaled streamwise, (b) wall-normal and (c) spanwise
velocity variances, and (d) outer-scaled Reynolds stresses profile. The legend
in (a) holds for all panels.

Figure 6 presents the second-order moments of the fluid and particle
velocities. Particle fluctuations are usually higher in the near-wall region
and smaller or similar to that of the fluid in the bulk. The modifications
with respect to the unladen case are much smaller for cases VD and SD,
compared with case SD. As discussed in Costa et al. (2020), the minor
deviations of VD with respect to the unladen case are attributed to a lower
statistical convergence of the massive particle-resolved simulation, and a
slight two-way coupling effect, while the more pronounced differences in case
D are caused by the increase in mean wall shear. Compared with these two
cases, case SD shows three significant differences: a strong reduction of the
streamwise velocity fluctuations near the wall and enhancement in the bulk, a
large decrease of the wall-normal and spanwise velocity variances across the
channel, and a reduction in Reynolds shear stresses. This observation is
consistent with particle-modelled DNS and experiments of particle-laden wall-
bounded turbulence (see, e.g. Kulick et al., 1994; Zhao et al., 2010; Richter,
2015; Capecelatro et al., 2018; Wang et al., 2019). Focusing on the fluid
Reynolds shear stress, case D features a slight increase of the peak value and
slope in the outer region, while the opposite applies to case SD. Note that
the particle Reynolds stresses are always larger than the fluid ones.

Despite the overall increase in drag, the significant decrease of the Reynolds
shear stresses for case SD is consistent with the reduced contribution to the
overall drag from the fluid turbulence shown in figure 2 by the friction
Reynolds number based on the velocity based on the slope of the Reynolds
stresses, $u_{\tau}^{t}=\sqrt{\partial_{y/h}\left\langle
u^{\prime}v^{\prime}\right\rangle|_{y/h=1}}$. This quantity initially
increases with the mass loading (case VD to D), denoting an enhancement of
turbulence and its induced drag. Further increasing the solid mass fraction,
(case SD), $Re_{\tau}^{t}$ decreases – even below the single-phase flow –
while the conventional friction Reynolds number increases. In other words,
while turbulence is attenuated, a fundamentally different mechanism acts to
anyway increase the drag. Also, despite showing drag increase, case D features
statistics that highly resemble those of the single-phase flow. Hence, only
case SD shows a more intricate two-way coupling mechanism, at play over the
whole channel region.

To gain further insight, we examine the streamwise momentum balance of the
fluid,

${\tau}={\rho}\left(1-\left\langle\phi\right\rangle\right)\left(\nu\frac{\mathrm{d}\left\langle
u\right\rangle}{\mathrm{d}y}-\left\langle
u_{f}^{\prime}v_{f}^{\prime}\right\rangle\right)+\left\langle\phi\right\rangle\tau_{p}=\rho
u_{\tau}^{2}\left(1-\frac{y}{h}\right)\mathrm{,}$ (2)

where the first and second terms denote the fluid viscous and Reynolds
stresses, while the last term the total particle stress
$\tau_{p}=\tau_{p,\nu}-\rho_{p}\left\langle
u_{p}^{\prime}v_{p}^{\prime}\right\rangle$, with $\tau_{p,\nu}$ including the
viscous and collisional particle stresses. The data are shown in figure
7(a-c), where the insets depict the local relative contribution of each term
to the total momentum transfer; panel (d) reports the share of each momentum
transfer mechanism to the overall drag, normalized by the overall single-phase
drag. These contributions are related to the stress profiles via the FIK
(Fukagata, Iwamoto & Kasagi) identity (Fukagata et al., 2002; Yu et al.,
2021), with the contribution of each stress mechanism $\tau_{i}$ given by the
following weighted integral: $\smallint_{0}^{1}3(1-y)\tau_{i}\,\mathrm{d}y$.
Expectedly, the stress budget pertaining to case VD resembles that of the
single-phase flow, with negligible contribution from the particles. Case D
displays a noticeable, though small, contribution of the particle stresses
which have a peak of about $5\%$ close to the wall and are of the order of the
viscous stresses in the bulk. Finally, the total particle stresses show a
significant relative contribution to the total in case SD, of about $20\%$.
Examining panel (d) of the same figure, we note that the sum of the viscous
and turbulent contributions corresponds to the single-phase drag in case VD,
with the particle stress showing an almost vanishing two-way effect. In case
D, the sum of the viscous and turbulent stress contributions is higher than
the overall single-phase drag, which is further increased by the particle
stress. This indicates that the overall drag is increased by two different
mechanisms: a direct one – the particles induce an extra stress, localized
near the wall where the average slip velocity is not negligible; and an
indirect effect – the solid phase enhances the Reynolds stress and hence the
turbulence. Notably, the sum of the viscous and turbulent stress contributions
is smaller than the overall single-phase drag in case SD; the particle stress,
however, counterbalances this reduction and the overall drag is still higher
than the corresponding unladen case. This confirms that we have turbulent drag
reduction and a strong contribution of the total particle stress, which
combine to a net drag increase.

(VD)(D)
(SD)

Figure 7: Budget of streamwise momentum as a function of the outer-scaled
wall-normal distance, for the different cases. The inset shows the same
budget, but normalized by the local value of the total stresses. The bottom-
right panel shows the relative contribution of each term in the stress budget
to the mean wall friction, normalized by that of the unladen case, i.e.,
$\smallint_{0}^{1}3(1-y)\tau_{i}\,\mathrm{d}y/(\rho u_{\tau}^{2\,sph})$
(Fukagata et al., 2002; Yu et al., 2021), with $\tau_{i}$ the different
contributions, and $u_{\tau}^{2\,sph}$ the wall friction velocity of the
unladen case.

Let us take the limit of vanishing volume fraction, $\phi\to 0$, but finite
mass fraction $\psi=\rho_{p}\phi$, of equation eq. (2),

$\tau_{\mathrm{two-way}}=\rho\left(\nu\frac{\mathrm{d}\left\langle
u\right\rangle}{\mathrm{d}y}-\left\langle
u_{f}^{\prime}v_{f}^{\prime}\right\rangle\right)-\left\langle\psi\right\rangle\left\langle
u_{p}^{\prime}v_{p}^{\prime}\right\rangle\approx\rho
u_{\tau}^{2}\left(1-\frac{y}{h}\right)\mathrm{.}$ (3)

This limit corresponds to negligible particle excluded volume, but finite
effects of the particle mass, i.e. two-way coupling conditions. In this limit,
the contribution from the correlated particle fluctuations corresponds to a
particle direct contribution to the drag. When plotting the terms in equations
(3), see figure 8, we note that the difference between the total stress and
the sum of the different terms is very small, confirming that the contribution
of the particle inertial shear stress to the budget is dynamically significant
and increases with the mass fraction. Moreover, the “two-way coupling” budget
explains the true nature of the turbulence attenuation for case SD. As the
mass loading is increased, the fluid Reynolds stresses progressively give in
to particle Reynolds stresses, which alter the nature of the flow. A
significant portion of the total stress is spent to accelerate the inertial
particles, which in turn lowers the fluid Reynolds stresses. Consequently, the
flow shows fluid turbulence attenuation. Unlike correlated fluid–fluid or
fluid–particle velocity fluctuations, particle–particle correlated motions do
not sustain near-wall turbulence, as demonstrated in Capecelatro et al. (2018)
for the fluid turbulence kinetic energy budget.

(D)(SD)

Figure 8: Budget of streamwise momentum in the two-way coupling limit, for
cases D and SD, (see equation (3)). The dashed grey line corresponds to the
difference between the sum of the terms in equation (3), and the total stress,
to assess the validity of the approximation. Case VD was omitted from this
figure, since given the negligible two-way coupling effect the budget is
trivial, resembling that of the single-phase flow.

Since the net effect of the particle Reynolds stresses is an increased drag
for case SD, we expect that at higher mass fractions the drag may eventually
decrease – the contribution of Reynolds stresses to the momentum transport
would decrease and the flow becomes smoother and eventually laminarizes.
Consequently, the correlated particle velocity fluctuations would also
decrease, drastically reducing the overall flow drag. We have actually
observed such a laminarizing and drag-reducing trend in a preliminary
simulation of a case with the same governing parameters as SD, but $10$ times
the number of particles (note that here particle excluded volume effects may
grow important). Note that, reductions of the skin friction coefficient have
also been observed in the two-way coupling point-particle DNS by Vreman
(2007); Zhao et al. (2010) and in the volume-filtered Euler–Lagrange
simulations of Capecelatro et al. (2018), with the latter study observing
laminarization at sufficiently high mass fraction. In particular, Capecelatro
et al. (2018) reported that the fluid velocity fluctuations were re-energized
at even higher mass fractions, due to the growing importance of turbulence-
producing particle–fluid velocity correlations. New experiments and fully
resolved simulations are hence needed to further support and explain these
observations at higher mass fractions.

##  4 Conclusions

We have used particle-resolved DNS to study the near-wall turbulence
modulation by small inertial particles. Three cases have been considered with
volume fractions progressively increased by one order of magnitude, chosen in
the one-way and two-way coupling regimes. The two densest cases show a non-
negligible turbulence modulation, however, of a fundamentally different
nature. The case with $\Psi=3.4\%$ ($\Phi=0.034\%$) features drag increase due
to increased inertia near the wall, with increased Reynolds stresses due to
particles travelling at high particle-to-fluid slip velocity, and dynamics
close to that of the single-phase flow at slightly higher Reynolds number.
Increasing the volume and mass fraction by a factor of $10$, the particle
presence is felt over the entire channel, however, the relative drag increase
is significantly reduced with respect to that of the previous configuration.
The streamwise momentum balance in the two-way coupling limit of vanishing
volume fraction, but finite mass fraction, reveals that correlated particle
fluctuations seize a larger share of the total stresses as the mass fraction
increases. This results in a reduction of the fluid Reynolds stresses, and
consequently of the turbulent drag. For the parameter setting considered here,
this turbulent drag reduction is counterbalanced by the contribution of
correlated particle velocity fluctuations; however, a further increase in mass
fraction may lead to a net drag reduction, something which has been noticed in
previous particle-modelled DNS, and should inspire future experimental and
particle-resolved numerical studies.

This work was supported by the European Research Council grant no.
ERC-2013-CoG-616186, TRITOS, the grant BIRD192032/19 from the University of
Padova, and the University of Iceland Recruitment Fund grant no. 1515-151341,
_TURBBLY_. We acknowledge the computing time provided by SNIC (Swedish
National Infrastructure for Computing), and PRACE for awarding us access to
the supercomputer Marconi, based in Italy at CINECA under project 2017174185 –
DILPART.

Declaration of Interests. The authors report no conflict of interest.

## References

  * Balachandar & Eaton (2010) Balachandar, S & Eaton, John K 2010 Turbulent dispersed multiphase flow. Annu. Rev. Fluid Mech. 42, 111–133.
  * Battista et al. (2019) Battista, F, Mollicone, J-P, Gualtieri, P, Messina, R & Casciola, CM 2019 Exact regularised point particle (erpp) method for particle-laden wall-bounded flows in the two-way coupling regime. J. Fluid Mech. 878, 420–444.
  * Breugem (2012) Breugem, W-P 2012 A second-order accurate immersed boundary method for fully resolved simulations of particle-laden flows. J. Comput. Phys. 231 (13), 4469–4498.
  * Capecelatro et al. (2018) Capecelatro, Jesse, Desjardins, Olivier & Fox, Rodney O. 2018 On the transition between turbulence regimes in particle-laden channel flows. J. Fluid Mech. 845, 499–519.
  * Costa (2018) Costa, P 2018 A fft-based finite-difference solver for massively-parallel direct numerical simulations of turbulent flows. Comput. Math. Appl. 76 (8), 1853 – 1862.
  * Costa et al. (2015) Costa, P, Boersma, BJ, Westerweel, J & Breugem, W-P 2015 Collision model for fully resolved simulations of flows laden with finite-size particles. Phys. Rev. E 92 (5), 053012.
  * Costa et al. (2020) Costa, P, Brandt, L & Picano, F 2020 Interface-resolved simulations of small inertial particles in turbulent channel flow. J. Fluid Mech. 883.
  * Crowe et al. (1977) Crowe, CT, Sharma, MPt & Stock, DE 1977 The particle-source-in cell (psi-cell) model for gas-droplet flows. J. Fluids Engng. 99 (2), 325–332.
  * Fornari et al. (2016) Fornari, W, Formenti, A, Picano, F & Brandt, L 2016 The effect of particle density in turbulent channel flow laden with finite size particles in semi-dilute conditions. Phys. Fluids 28 (3), 033301.
  * Fröhlich et al. (2018) Fröhlich, K, Schneiders, L, Meinke, M & Schröder, W 2018 Validation of lagrangian two-way coupled point-particle models in large-eddy simulations. Flow Turbul. Combust. 101 (2), 317–341.
  * Fukagata et al. (2002) Fukagata, Koji, Iwamoto, Kaoru & Kasagi, Nobuhide 2002 Contribution of reynolds stress distribution to the skin friction in wall-bounded flows. Phys. Fluids 14 (11), L73–L76.
  * Gualtieri et al. (2015) Gualtieri, P, Picano, F, Sardina, G & Casciola, CM 2015 Exact regularized point particle method for multiphase flows in the two-way coupling regime. J. Fluid Mech. 773, 520–561.
  * Horne & Mahesh (2019) Horne, WJ & Mahesh, K 2019 A massively-parallel, unstructured overset method to simulate moving bodies in turbulent flows. J. Comput. Phys. 397, 108790.
  * Horwitz & Mani (2016) Horwitz, JAK & Mani, A 2016 Accurate calculation of stokes drag for point–particle tracking in two-way coupled flows. J. Comput. Phys. 318, 85–109.
  * Ireland & Desjardins (2017) Ireland, PJ & Desjardins, O 2017 Improving particle drag predictions in euler–lagrange simulations with two-way coupling. J. Comput. Phys. 338, 405–430.
  * Kim & Moin (1985) Kim, J & Moin, P 1985 Application of a fractional-step method to incompressible navier-stokes equations. J. Comput. Phys. 59 (2), 308–323.
  * Kuerten & Vreman (2016) Kuerten, J.G.M. & Vreman, A.W. 2016 Collision frequency and radial distribution function in particle-laden turbulent channel flow. Int. J. Multiphase Flow 87, 66–79.
  * Kulick et al. (1994) Kulick, Jonathan D, Fessler, John R & Eaton, John K 1994 Particle response and turbulence modification in fully developed channel flow. J. Fluid Mech. 277, 109–134.
  * Marchioli & Soldati (2002) Marchioli, Cristian & Soldati, Alfredo 2002 Mechanisms for particle transfer and segregation in a turbulent boundary layer. J. Fluid Mech. 468, 283–315.
  * Mehrabadi et al. (2018) Mehrabadi, M, Horwitz, JAK, Subramaniam, S & Mani, A 2018 A direct comparison of particle-resolved and point-particle methods in decaying turbulence. J. Fluid Mech. 850, 336–369.
  * Pakseresht et al. (2020) Pakseresht, Pedram, Esmaily, Mahdi & Apte, Sourabh V. 2020 A correction scheme for wall-bounded two-way coupled point-particle simulations. J. Comput. Phys. 420, 109711\.
  * Picano et al. (2015) Picano, F, Breugem, W-P & Brandt, L 2015 Turbulent channel flow of dense suspensions of neutrally buoyant spheres. J. Fluid Mech. 764, 463–487.
  * Richter (2015) Richter, David H 2015 Turbulence modification by inertial particles and its influence on the spectral energy budget in planar couette flow. Phys. Fluids 27 (6), 063304.
  * Rouson & Eaton (2001) Rouson, DWI & Eaton, JK 2001 On the preferential concentration of solid particles in turbulent channel flow. J. of Fluid Mech. 428, 149.
  * Sardina et al. (2012) Sardina, G, Schlatter, P, Brandt, L, Picano, F & Casciola, CM 2012 Wall accumulation and spatial localization in particle-laden wall flows. J. Fluid Mech. 699 (1), 50–78.
  * Schneiders et al. (2017) Schneiders, L, Meinke, M & Schröder, W 2017 Direct particle–fluid simulation of kolmogorov-length-scale size particles in decaying isotropic turbulence. J. Fluid Mech. 819, 188–227.
  * Sundaram & Collins (1997) Sundaram, S & Collins, LR 1997 Collision statistics in an isotropic particle-laden turbulent suspension. part 1. direct numerical simulations. J. Fluid Mech. 335, 75–109.
  * Tiederman et al. (1985) Tiederman, WG, Luchik, Thomas S & Bogard, DG 1985 Wall-layer structure and drag reduction. J. Fluid Mech. 156, 419–437.
  * Uhlmann (2005) Uhlmann, M 2005 An immersed boundary method with direct forcing for the simulation of particulate flows. J. Comput. Phys. 209 (2), 448–476.
  * Vreman (2007) Vreman, A. W. 2007 Turbulence characteristics of particle-laden pipe flow. J. Fluid Mech. 584, 235–279.
  * Wang et al. (2019) Wang, G, Fong, KO, Coletti, F, Capecelatro, J & Richter, DH 2019 Inertial particle velocity and distribution in vertical turbulent channel flow: A numerical and experimental comparison. Int. J. Multiphase Flow 120, 103105.
  * Yu et al. (2021) Yu, Zhaosheng, Xia, Yan, Guo, Yu & Lin, Jianzhong 2021 Modulation of turbulence intensity by heavy finite-size particles in upward channel flow. J. Fluid Mecha. 913.
  * Zhao et al. (2010) Zhao, LH, Andersson, H I & Gillissen, JJJ 2010 Turbulence modulation and drag reduction by spherical particles. Phys. Fluids 22 (8), 081702.