The exact regularized point particle method is used to characterize the turbulence modulation in two-way momentum-coupled direct numerical simulations of a turbulent pipe flow. The turbulence modification is parametrized by the particle Stokes number, the mass loading, and the particle-to-fluid density ratio. The data show that in the wide region of the parameter space addressed in the present paper, the overall friction drag is either increased or unaltered by the particles with respect to the uncoupled case. In the cases where the wall friction is enhanced, the fluid velocity fluctuations show a substantial modification in the viscous sub-layer and in the buffer layer. These effects are associated with a modified turbulent momentum flux toward the wall. The particles suppress the turbulent fluctuations in the buffer region and concurrently provide extra stress in the viscous sub-layer. The sum of the turbulent stress and the extra stress is larger than the Newtonian turbulent stress, thus explaining the drag increase. The non-trivial turbulence/particles interaction turns out in a clear alteration of the near-wall flow structures. The streamwise velocity streaks lose their spatial coherence when two-way coupling effects are predominant. This is associated with a shift of the streamwise vortices toward the center of the pipe and with the concurrent presence of small-scale and relatively more intense vortical structures near the wall.

## I. INTRODUCTION

Particle-laden turbulent flows are ubiquitous in many technological applications and natural phenomena. Among them there are the transport of sediments in rivers,^{1} the dynamics of ashes and plumes in the atmosphere,^{2} the particle turbulent transport, and their deposition and entrainment in turbulent boundary layer,^{3} chemical plants' applications,^{4} or the genesis and evolution of small droplets.^{5}

Given the complex interaction between the carrier and the disperse phase, a multi-scale and multi-physics approach is necessary to model the particle/turbulence interaction.^{6–9} The presence of flow confinement introduces additional complexity. The preferential segregation of the particles near the walls^{10–12} leads the particles to affect the turbulence in the buffer layer.^{13} It follows a modification of the regeneration cycle of the vortical structures and of the flow topology,^{14–17} which ultimately manifests, at least for finite-size neutrally buoyant particles, in a modification of the well-known law of the wall.^{18}

In the context of Eulerian–Lagrangian approaches—for Eulerian–Eulerian simulations the reader is referred to Ref. 19—multiphase flows are simulated with two main approaches:^{20} (i) the resolved particle simulations that resolve the particle length-scale, e.g., diameter, on the computational grid,^{21–26} and (ii) the point-particle simulations that model the particles as material points. The former approach, even though is based on first principles, cannot be pursued when the particles are relatively small with respect to the smallest hydrodynamical length-scale, or when the particle number is large, due to its tremendous computational cost. However, when the particle diameter *d _{p}* is comparable with the smallest hydrodynamical length-scale $ \u2113 *$, it is reasonable to adopt the point-particle method.

^{27}Given the small value of the particle Reynolds number, the hydrodynamic force can be parametrized in terms of the Stokes drag.

^{28,29}In conditions where the particle volume fraction $ \varphi V$ (ratio between the particle volume and the fluid volume) is small, the fluid simply carries the disperse phase and it is not affected by the presence of the particles (one-way coupling regime). However, when the particle-to-fluid density ratio $ \rho p / \rho f$ is large, there exist conditions where the volume fraction is still small to neglect inter-particles collisions, but the mass loading $ \varphi = ( \rho p / \rho f ) \u2009 \varphi V$ is not negligible.

^{30}In these conditions, the disperse phase affects the momentum dynamics of the carrier phase (two-way coupling regime).

Within the point-particle approximation, the force of the particles on the fluid is concentrated at the actual particle position by a Dirac delta function and needs to be regularized. The simplest regularization is provided by the particle-source in cell (PIC) method.^{27} The approach, however, has several drawbacks, mainly: (i) the solution is not grid convergent depending on how many particles per computational cell are available and (ii) it is not possible to evaluate the unperturbed fluid velocity at the particle position (needed to evaluate the Stokes drag) since every single particle locally modifies the fluid velocity.^{31,32} The results based on the PIC approach are quite scattered. In turbulent channel flows, both drag reduction and turbulent fluctuation attenuation are observed^{33} with $ \rho p / \rho f = 1024 , \u2009 S t + = 30$, and $ \varphi = 0.96$, see also Ref. 34 where the mass loading is changed from 0 to 0.96. Other simulations^{35} show an overall increase in the turbulent fluctuations for relatively small particles and a decrease for large particles; here, the particle Stokes number $ S t +$ spans the range $ 0.5 \u2013 125$ and the density ratio varies from 35 to 8650 at fixed mass loading $ \varphi = 0.3$. Variable outcomes are obtained depending on the mass loading of the suspension and on the particle Stokes number in Ref. 36 ( $ \varphi = 0.2 , 0.4 , 2 , \u2009 S t + = 60 , 190$, and $ \rho p / \rho f = 2083 , 7333$). Turbulent pipe flow simulations^{37,38} with $ \varphi = 0.22 , 0.44 , 0.89 , \u2009 S t + = 14 , 29 , 58$, and $ \rho p / \rho f = 500 , 1000 , 2000$ provide similar results as those observed in the channel flow, while in turbulent boundary layers, an increase in the skin friction coefficient is observed.^{39}

In the context of the Euler–Lagrangian approaches, many techniques have been developed from the numerical point of view to get rid of the particle self-disturbance.^{40–47} In the exact regularized point particle method (ERPP),^{46} the inter-phase momentum coupling occurs by a mechanism of vorticity generation by the particles and its diffusion by viscosity. This mechanism allows a physically consistent description of the disturbances produced by the particles in the framework of the unsteady Stokes flows. Moreover, the particle self-disturbance is known and can be removed from the background fluid velocity when the Stokes drag is evaluated. The approach, which provides convergent turbulent statistics,^{48,49} has been generalized to deal with wall bounded flows.^{50} Recently, it has also been used in different contexts, e.g., polymer-laden flows,^{51} and to model several multi-physical phenomena, e.g., evaporating droplets in shear flows,^{52} in isotropic regions of turbulent sprays,^{53} and in reacting flows.^{54} Concerning evaporating droplets, Ref. 55 would be a valuable benchmark to compare the point-particle approach with a finite-size droplet simulation.

Experimentally, the motion, deposition, entrainment, and spatial distribution of the particles in boundary layers and channel flows have been addressed in Refs. 56 and 57. The study reports larger velocity gradients close to the wall, which corresponds to a larger wall shear stress with respect to the unladen case. The augmentation of the turbulent velocity fluctuations close to the wall is also observed. Other experimental measurements show an increase in friction and in turbulent fluctuations at the wall.^{58–61} In the pipe flow, an increase in the wall shear stress and in the turbulent fluctuations near the wall is also observed.^{62–65}

The scattering in the results is unavoidably due to the complex multi-scale interaction between the carrier and the disperse phase, which is controlled by many parameters, e.g., the turbulence Reynolds number, the particle Stokes number, the particle-to-fluid density ratio, the mass loading (or volume fraction) of the suspension, and the particle Reynolds number. In any case, both resolved particle direct numerical simulations and the experimental results seem to agree on the fact that the addition of the particles increases the wall shear stress with respect to the unladen case (drag increase) and that the turbulent fluctuations are augmented in the near-wall region.

The present paper investigates the modification of wall turbulence in a turbulent pipe flow laden with inertial particles. The analysis aims to explore a wide range of the phase space parametrized by the particle Stokes number, the mass-loading, and the particle-to-fluid density ratio, at fixed turbulence Reynolds number. This work is a continuation/extension of a previous study^{50} with a deeper and more extensive analysis, of both fluid and particle velocity fluctuations statistics. The impact of the particles' feedback on the flow structures of wall-bounded flows, i.e., the velocity streaks and the coherent quasi-streamwise vortices, is also addressed. The inter-phase momentum coupling is exploited by the generalized ERPP approach.^{46,50} The issue of whether the particle feedback produces an overall increase or decrease in the friction drag and whether the turbulent fluctuations are augmented or reduced by the particles is discussed. In the region of the parameter space covered by the present simulations, the data show that the particles always increase or at most leave unaltered the friction drag and that the fluid velocity fluctuations are augmented in the viscous sub-layer. The statistics of the dispersed phase in terms of velocity fluctuations follow the same behavior as the fluid phase. The non-trivial turbulence/particle interaction turns out in a clear alteration of the streamwise velocity streaks that lose their spatial coherence due to two-way coupling effects. Concurrently, the streamwise vortices are shifted toward the center of the pipe, and new small-scale and relatively intense vortical structures appear in the near-wall region.

The paper is organized as follows: Sec. II briefly summarizes the generalized ERPP approach and reports the simulation setup and parameters for the turbulent pipe flow. Section II discusses the turbulence modification and the particle statistics in the two-way coupling regime. Section II C presents the conclusions of this work.

## II. PARTICLE LADEN TURBULENT PIPE FLOW

*θ*, radial

*r*, and axial

*z*dimensions are normalized with the pipe radius

*R*. Periodic boundary conditions are applied in the axial and azimuthal direction, and impermeability and no-slip conditions are enforced at the walls. In Eq. (1), $ Re b , 0 = \rho f U b , 0 R / \mu $ is the bulk Reynolds number, $ U b , 0 = Q 0 / ( \pi R 2 )$ is the bulk velocity, and

*Q*

_{0}is the flow rate of the reference uncoupled case (no particle back-reaction). As usual,

*μ*is the dynamic viscosity,

*ρ*is the fluid density, and $ \nu = \mu / \rho f$ is the kinematic viscosity. The flow is sustained by a constant mean pressure gradient that is applied in the direction of the axial unit vector $ e z$. The dimensionless pressure is expressed as $ P = d p / d z | 0 ( z \u2212 z 0 ) + p ( \theta , r , z , t )$. Equations (1) are solved in cylindrical coordinates exploiting a second-order finite difference discretization on a staggered grid. The classical Chorin's projection method

_{f}^{66}is used to enforce the divergence-free constraint imposed by the mass balance. Both convective and diffusive terms are explicitly integrated in time using a third-order low-storage Runge–Kutta method. In the context of the exact regularized point particle method, the inter-phase momentum coupling is given by the field $ f ( x , t )$ in Eq. (1),

*ϵ*that regularizes the particles' feedback. The tilde denotes the particle image to account for the presence of the wall, and $ D p$ denotes the hydrodynamic force, i.e., the Stokes drag. The reader is referred to the original papers

^{46,50}for the complete, thorough derivation of the feedback term.

*τ*is the particle relaxation time. When the particle Reynolds number is small and the particle-to-fluid density ratio is large, the hydrodynamic force reduces to the Stokes drag and the related Faxen correction.

_{p}^{28,29}A purely elastic collision model is used to manage the particle/wall interaction. In the expression of the hydrodynamic force, the hat denotes the fluid velocity evaluated at the particle position in the absence of the $ p th$ particle, i.e., it is the

*undisturbed*velocity field. This field is evaluated by subtracting from the background flow field the

*particle self-disturbance*that is known in a closed form in the ERPP approach. Equations (3) are integrated in time with the same Runge–Kutta method employed for the carrier phase.

In wall turbulence, it is usual to consider inner or wall units, i.e., the friction velocity $ u * = \tau w / \rho f$, where *τ _{w}* is the average wall shear stress, the viscous length $ \u2113 * = \nu / u *$, and the viscous timescale $ \tau * = \u2113 * / u *$. The distance from the pipe wall in inner units is denoted by $ y + = ( 1 \u2212 r / R ) Re *$, where $ Re * = u * R / \nu $ is the friction Reynolds number. The same distance in external units is denoted by $ y = 1 \u2212 r / R$. All the simulations are performed with the same friction Reynolds number $ R e * = 180$, corresponding to a bulk Reynolds number of $ R e b , 0 = 2650$ in the uncoupled case. The grid resolution is $ N \theta \xd7 N r \xd7 N z = 576 \xd7 129 \xd7 576$ in the azimuthal, wall-normal, and axial directions, respectively. In the radial direction, the grid is clustered near the wall with a minimum spacing of $ \Delta r + | w = 0.5$ which gradually increases toward the centerline reaching $ \Delta r + | 0 = 2$. The grid resolution in the azimuthal and axial directions is $ ( R \Delta \theta ) + = 3$ and $ \Delta z + = 3$, respectively.

The disperse phase is characterized by the inner-scale Stokes number defined as $ St + = \tau p / \tau * = St b Re * 2 / Re b , 0$. The mass loading of the suspension is the ratio between the total mass of the disperse and carrier phases, that is, $ \varphi = ( \rho p / \rho f ) \varphi V$, where $ \varphi V = N p \u2009 V p / V f$ is the volume fraction, where *V _{p}* is the volume of the single particle, and

*V*is the volume of the fluid domain $D$.

_{f}In summary, the dynamics of the suspension is controlled by a set of four dimensionless parameters $ { R e * ; \u2009 S t + ; \u2009 \varphi ; \u2009 \rho p / \rho f}$ that correspond to the following physical assumptions: (i) the volume fraction is small (dilute suspension) to neglect particle–particle interactions, (ii) the particle diameter $ d p +$ is of the order of the viscous length and the particle Reynolds number is small, and (iii) the density ratio $ \rho p / \rho f$ is sufficiently large to allow non-negligible mass loading at small volume fractions. The parameters of the different cases are summarized in Table I. The simulations are grouped into three sets. In the first set, the mass loading $\varphi $ is changed keeping the Stokes number and density ratio fixed. The second set addresses the effects of the Stokes number at fixed mass loading and density ratio. Finally, the density ratio is changed at fixed mass loading and Stokes number.

$\varphi $ . | $ \rho p / \rho f$ . | $ S t +$ . | St
. _{b} | $ d p +$ . | N
. _{p} |
---|---|---|---|---|---|

0 | ⋯ | ⋯ | ⋯ | ⋯ | ⋯ |

0.2 | 180 | 10 | 0.82 | 1 | 244 290 |

0.4 | 180 | 10 | 0.82 | 1 | 488 580 |

0.6 | 180 | 10 | 0.82 | 1 | 732 870 |

0.4 | 180 | 15 | 1.23 | 1.23 | 265 950 |

0.4 | 180 | 20 | 1.64 | 1.41 | 172 739 |

0.4 | 180 | 80 | 6.54 | 2.82 | 21 592 |

0.4 | 90 | 10 | 0.82 | 1.41 | 345 479 |

0.4 | 360 | 10 | 0.82 | 0.70 | 690 957 |

0.4 | 560 | 10 | 0.82 | 0.57 | 861 775 |

$\varphi $ . | $ \rho p / \rho f$ . | $ S t +$ . | St
. _{b} | $ d p +$ . | N
. _{p} |
---|---|---|---|---|---|

0 | ⋯ | ⋯ | ⋯ | ⋯ | ⋯ |

0.2 | 180 | 10 | 0.82 | 1 | 244 290 |

0.4 | 180 | 10 | 0.82 | 1 | 488 580 |

0.6 | 180 | 10 | 0.82 | 1 | 732 870 |

0.4 | 180 | 15 | 1.23 | 1.23 | 265 950 |

0.4 | 180 | 20 | 1.64 | 1.41 | 172 739 |

0.4 | 180 | 80 | 6.54 | 2.82 | 21 592 |

0.4 | 90 | 10 | 0.82 | 1.41 | 345 479 |

0.4 | 360 | 10 | 0.82 | 0.70 | 690 957 |

0.4 | 560 | 10 | 0.82 | 0.57 | 861 775 |

## III. RESULTS

### A. Fluid statistics: Mean velocity and velocity fluctuations profiles

Figure 1 shows the mean velocity profiles in semi-logarithmic scale in internal units (left column) and in linear scale in external units (right column). Panels (a) and (b) show the effect of the mass loading, panels (c) and (d) the effect of the Stokes number, and panels (e) and (f) the effect of the density ratio. The drag increase manifests as a decrease in the flow speed since the pressure gradient is the same in all the studied cases. When rescaled in internal units, the velocity profile in the viscous sub-layer is not modified by the inter-phase momentum coupling, while the buffer and logarithmic regions are significantly altered, especially when a notable drag increase occurs. The present results, obtained with a point-particle approach, are capturing the important physical behavior of a particle-laden flow, that is, the increase in the overall drag. Previously, similar results have only been achieved via resolved particle simulations^{23,67,68} or by exploiting the volume-filtered Navier–Stokes approach to account for excluded volume effects.^{42,69} Similar behavior also occurs at low (order one) particle-to-fluid density ratio, as documented both numerically^{22,24} and experimentally.^{58–60,62} The present results document that the increase in the drag is also present at higher values of the particle-to-fluid density ratios and at small values of the Stokes number.

The increase in the friction drag manifests as a depletion of the mass flow rate. Figure 2 shows the bulk velocity for all the simulations and summarizes at a glance the effect that the different parameters have on the bulk velocity modification. The increase in the drag is observed when the mass loading is increased, when the Stokes number is relatively small, or when the density ratio is order 100. This latter circumstance is not trivial and roughly corresponds to actual cases such as, for example, when medicinal powders are used for inhalable drug delivery, when carbon dust is transported, in the food industry when powders result from the processing of cereals, and when sawdust is produced in wood manufacturing. Additionally, when the carrier fluid is water or a relatively dense fluid, common materials turn out to have a relatively small density ratio. Note that *U _{b}* is a function of $ [ R e * , S t + , \varphi , \rho p / \rho f ]$ and Fig. 2 provides only three one-dimensional views at fixed $ R e *$ of this four-variable phase space.

The characterization of the turbulence modulation is completed by addressing the mean profiles of velocity variances and Reynolds shear stresses. The axial and radial mean velocity variance plotted against the wall-normal distance is shown in Figs. 3 and 4, respectively. The largest effects on the flow are seen in the viscous sub-layer and in the buffer layer. In the bulk, a substantial modification is also observed for the radial (wall normal) velocity component. In the buffer region, two-way coupling effects decrease the peak of velocity variances and lead to an augmentation of the variances in the viscous sub-layer ( $ 2 < y + < 5$). This scenario corresponds to cases where a substantial drag increase is observed. The effect of the mass loading is significant since the peak in the buffer region is attenuated and a new peak appears in the viscous sub-layer. The departure from the classical Newtonian behavior is attenuated when the Stokes number or the density ratio is increased. At larger Stokes numbers and/or density ratios, the viscous sub-layer is unaffected even though the depletion of the peak in the buffer region is still present.

Figure 5 completes the analysis of velocity variances by showing the Reynolds shear stress. By increasing the mass loading, the peak of the Reynolds stresses in the buffer region is attenuated, while an increase appears in the viscous sub-layer. The same occurs at a relatively small Stokes number, while when the Stokes number is progressively increased, the shear stress profile recovers the behavior of the one-way coupled simulation except for the intensity of the peak in the buffer region, which is still below the uncoupled case. A similar effect is observed by increasing the density ratio.

^{70}The result is a balance equation that relates the pressure drop to the flow rate via the turbulent stress and the particles' extra stress, which takes part in the budget because of two-way coupling effects. The budget reads

*U*. In the turbulent case, the presence of the Reynolds shear stress results in a reduced flow rate for a given pressure drop. In the two-way coupling regime, the extra stress also plays a role in determining the flow rate—together with the (modified) turbulent stress—since it can absorb part of the available pressure drop. Note that both the turbulent stress,

_{b}*τ*, and the extra-stress,

_{t}*τ*, enter in the budget through the weight factor

_{e}*r*

^{2}. It follows that any slight modification of the stresses in the near-wall region has a striking effect on the flow rate. The budget (4) suggests the analysis of the weighted stress profiles, $ r 2 \tau t$ and $ r 2 \tau e$, see Fig. 6. In the case at $ \varphi = 0.4 , \u2009 S t + = 10$, and $ \rho p / \rho f = 180$, panel (a), the augmentation of the Reynolds shear stress in the near-wall region and the presence of the extra stress explain the change in the profile of the cumulative turbulent stress $ r 2 ( \tau t + \tau e )$, which overall follows almost everywhere the data of the one-way coupling regime $ r 2 \tau t | \varphi = 0$, except in the near-wall region where it is strongly augmented. In contrast, for the particle population at $ \varphi = 0.4 , \u2009 S t + = 80$, and $ \rho p / \rho f = 180$ shown in panel (b), both the profiles of $ r 2 \tau t$ and $ r 2 \tau e$ are modified by the backreaction but their sum meets the profile $ r 2 \tau t | \varphi = 0$ of the uncoupled case. This explains why the drag is not modified even though particles carry significant extra stress. Finally, in the case at $ \varphi = 0.4 , \u2009 S t + = 10$, and $ \rho p / \rho f = 560$, panel (c), the absence of drag increase is related to the almost absent extra stress and a negligible modification of the turbulence Reynolds stress. In general, a lower peak of the Reynolds shear stresses in the buffer layer would result in drag reduction. However, the increase in Reynolds stress in the viscous sub-layer combined with the particle extra stress close to the wall largely overwhelms this effect, since the contribution of the viscous sub-layer, see the weighting factor

*r*

^{2}, is dominant. This leads to an overall increase in drag, associated with a higher compound turbulent stress.

A similar behavior, obtained here by modeling particles as concentrated masses, reproduces a turbulence regime that has been observed only in resolved particle simulations, see, e.g., Ref. 68. This indicates that the drag increase is not related to the particle's finite size but to the feedback itself, as long as it is modeled in a physically consistent way.

### B. Particle statistics: Mean velocity and velocity variances

The effect of the particle feedback on wall turbulence is expected to be larger where the particles tend to accumulate in space. The phenomenon of preferential segregation near a solid wall is known as turbophoresis and is clearly observed in the one-way coupled regime.^{71–73} In the two-way coupling regime, the phenomena still occur as discussed in Ref. 50. In the following, we extend the analysis to the particle mean velocity and velocity fluctuations.

Figure 7 shows the particle mean velocity compared to the fluid mean velocity. In panel a), the mass loading is changed, while the other parameters are fixed. At low mass loading, two-way coupling effects are small and both the fluid and the particle mean velocity match the behavior of the one-way coupled simulations. When $\varphi $ is increased, both fluid and particle velocity are substantially modified. Note that the particles share the same fluid mean velocity away from the wall, while when close to the buffer layer, the particle velocity lags behind the fluid velocity, indicating that the fluid is faster than the particles on average. Panel (b) addresses the effect of the Stokes number. The population with the higher Stokes number addressed, $ S t + = 80$, shows that the particle velocity lags behind the fluid velocity almost everywhere. Finally, panel (c) shows the effect of the density ratio where similar observations hold.

Figure 8 reports the particles' (symbols) and fluid (only for the one-way coupled case, solid line) velocity fluctuations. The plots in the left column address the axial velocity component, while the plots in the right column address the radial (wall-normal) velocity component. Panels (a) and (b) demonstrate the effect of mass loading. As $\varphi $ increases, the peak in the buffer region is progressively attenuated and velocity fluctuations correspondingly increase in the viscous-sub layer. This behavior is generic for both velocity components. Panels (c) and (d) demonstrate the effect of the Stokes number at fixed mass loading and density ratio. Concerning the axial particle velocity variance, the population with $ S t + = 10$ shows relatively intense fluctuations in the viscous sub-layer and a relatively small peak in the buffer region. When $ S t +$ increases, the particle velocity variance activates throughout the flow domain, while the radial velocity fluctuations show depletion. Finally, panels (e) and (f) address the effect of the density ratio and show how relatively light particle populations have large variances near the wall. In contrast, at higher particle-to-fluid density ratios, variances are particularly intense in the buffer region. The peak of the radial velocity fluctuations can be associated with particles that travel toward/away from the walls. Moreover, since the particles' and fluid mean velocity profiles are similar, the axial velocity fluctuations can be interpreted as a proxy of the momentum exchange between the two phases. The particles' shear stresses better illustrate this point.

Figure 9 shows the particle streamwise/wall-normal cross correlation, referred to as the particle Reynolds shear stresses. In the three panels of the figure, the Reynolds shear stresses of the fluid in the one-way coupling regime are also shown for comparison. Panel (a) addresses the effect of the mass loading and shows that at relatively small values of $\varphi $, the particle Reynolds shear stresses are particularly active in the buffer region. At increasing mass loading, the peak of the stresses is progressively eroded where particles' shear stress becomes more active in the viscous sub-layer. The effect of the Stokes number is addressed in panel (b) which shows that increasing $ S t +$, the stresses are progressively attenuated in the viscous sub-layer, while a definite peak appears in the buffer region. At the higher Stokes number, the peak is slightly shifted toward the outer part of the buffer layer and the overall intensity is depleted. Panel (c) presents the effect of the particle-to-fluid density ratio and shows that at increasing $ \rho p / \rho f$, the peak of the particles' shear stresses in the buffer region is increasingly more pronounced and the curve always lays above the corresponding fluid Reynolds shear stresses. The particle shear stress is related to a simple mechanism that involves the particles that come from the bulk of the flow and approach the walls (or vice versa). The shear is positive since the particle radial velocity fluctuations are positive (particles approach the wall), and they generate a positive streamwise velocity fluctuation since they bring a higher streamwise momentum. Close to the wall, the particle and the fluid mean velocity profile are similar, and it follows that the effect of this dynamics is to transfer positive streamwise momentum to the fluid (the particles push the fluid), hence generating high-velocity gradients at the wall, i.e., drag increase. Note that also for a particle traveling away from the wall, the shear stress is positive and results in a momentum exchange that goes from the fluid to the particles (the fluid pushes the particles). This generates a loop where the particle absorbs streamwise momentum from the fluid in the bulk and releases it near the wall, giving an alternative way of momentum transfer that ultimately increases the drag.

### C. Instantaneous near-wall structures

Turbulence and particle statistics are mostly affected with respect to the one-way coupling regime in the buffer layer and viscous sub-layer. This suggests analyzing the typical turbulence structures of this flow region, i.e., the streaks visualized by the instantaneous axial velocity fluctuations $ u \u2032 z$. The streaks are shown in Fig. 10, and the instantaneous particle configuration is superimposed below. Among the many cases addressed in this paper, only the cases corresponding to extreme values of the parameters of Fig. 2 are shown. The one-way coupled visualization is reported in panel (a) for comparison. When the particles' effect is negligible, that is at a high Stokes number, $ S t + = 80$ in panel (c) and high density ratio, $ \rho p / \rho f = 560$ in panel (e), the streaks show the same qualitative structure of the classical wall turbulence. Indeed, in these cases, the particles tend to segregate into high-speed streaks. However, when the momentum coupling is significant, that is for case $ \varphi = 0.6$ in panel (b), $ S t + = 10$ in panel (d) and $ \rho p / \rho f = 90$ in panel (f), the streaks are less evident. Their streamwise coherence appears to be lost, and the flow is characterized by large regions of negative velocity fluctuations and very confined regions of positive intense velocity fluctuations, see the red spots in the panels. The spots of positive velocity fluctuations correspond to the localized particle agglomerates whose concentration correlates at a glance with the velocity fluctuations peaks. Outside these regions, the particle distribution appears more uniform with respect to the other cases even though the instantaneous snapshots suggest that the particles' effect is more intense in the regions where clusters of particles form. The streaks are typically associated with the presence of coherent vortical structures. Figure 11 shows the iso-contours of the well-established *Q*-criterion, where $ Q = 0.5 ( | | \Omega | | 2 \u2212 | | S | | 2 ) , \u2009 \Omega = 0.5 ( \u2207 u \u2212 \u2207 u T )$, and $ S = 0.5 ( \u2207 u + \u2207 u T )$. A positive value of *Q* corresponds to regions of the flow field in which vorticity dominates the strain, i.e., the coherent vortical structures. The one-way coupling case is shown in panel (a) for comparison. Panel (b) corresponds to the case at $ S t + = 80$ where no substantial alteration of the drag is measured. Panels (c) and (d) report the structures for the case at $ \varphi = 0.6$, where, instead, the effect of the particle on the drag is significant. In all the panels, the isosurfaces of the *Q*-criterion are colored with the wall-normal distance. The case at *St* = 80 shows vortical structures qualitatively similar to the one-way coupling simulation, for the same value of the isosurface *Q* = 2. The same isosurface *Q* = 2 at $ \varphi = 0.6$ [panel (c)] shows that the coherent structures are moved in the center of the pipe and that a lawn of relatively small structures appears near the wall. These are highlighted by addressing the isosurface *Q* = 20, see panel (d), and surprisingly their wall normal distance matches the peak position in fluid velocity fluctuations. These small structures are very intense and are associated with the loss of spatial order of the near-wall streaks, which has been previously discussed.

## IV. FINAL REMARKS

The wall turbulence modulation due to small inertial particles has been investigated in a fully developed pipe flow exploiting the exact regularized point particle (ERPP) method to model the inter-phase momentum coupling.

Even though the results have been obtained with the drastic simplification of particles taken as point masses, it appears that in specific regions of the parameter space, the overall friction drag is increased by two-way coupling effects. In these cases, the fluid and the particle velocities show high-intensity fluctuations in the viscous sub-layer and an attenuation in the buffer region. Friction drag and fluctuations increase for the particle populations at relatively small Stokes numbers $ S t + = 10 \u2013 20$ and for the mass loading of $ \varphi = 0.4 \u2013 0.6$. The carrier and the disperse phase alteration is attenuated as the Stokes number increases beyond $ S t + = 80$. The particle-to-fluid density ratio also plays an important role. The turbulence modification is explained in terms of the alteration of the turbulent stresses. The particles provide extra stress in the viscous sub-layer. The sum of the turbulent stress and the extra stress is larger than the Newtonian turbulent stress, thus explaining the drag increase in terms of an alteration of the streamwise momentum flux toward the wall. The turbulence alteration is visualized at a glance by the modification of the near-wall streamwise velocity streaks that appear somehow destroyed in their streamwise coherence in cases where a significant drag increase is observed. This is associated with a shift of the streamwise vortices toward the center of the pipe and with the concurrent presence of small-scale and relatively more intense vortical structures near the wall.

## ACKNOWLEDGMENTS

The authors acknowledge CINECA for awarding access to supercomputing resource MARCONI based in Bologna, Italy, through ISCRA-B Project Nos. HP10BSK3XZ and HP10B8M78C, funded by the European Union. This work has been supported by Italian PNRR funds, CN-1 Spoke 6.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Francesco Battista:** Conceptualization (equal); Data curation (equal); Methodology (equal); Software (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). **Paolo Gualtieri:** Conceptualization (equal); Data curation (equal); Methodology (equal); Software (equal); Supervision (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). **Jean-Paul Mollicone:** Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). **Francesco Salvadore:** Software (equal); Writing – review & editing (equal). **Carlo Massimo Casciola:** Conceptualization (equal); Supervision (equal); Writing – review & editing (equal).

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.

## REFERENCES

*o*(

*N*) method for the numerical simulation of disperse systems: Potential flow of spheres

*IUTAM Symposium on Computational Approaches to Multiphase Flow*