1 Introduction

Particle-laden turbulent flows are commonly encountered in natural and industrial settings, such as sediment transport, the paper industry, pipeline transport, and fluidized beds. It is important to understand particle-turbulence interaction mechanisms to improve the macroscopic models for the multi-phase flows and the design of the related device. The point-particle-approximation-based direct numerical simulations (DNSs) have provided extensive insight into particle-turbulence interactions^[1-3]. However, in principle, the point-particle model is suited to situations where the particle size is smaller than the turbulence Kolmogorov length scale and the particle volume fraction is low. In recent years, interface-resolved DNS methods have been developed for particles larger than the Kolmogorov length scale^[4]. The essential features of the interface-resolved methods are that the interfaces between the particles and the fluid are resolved and the hydrodynamic forces on the particles are determined from the solution of the flow fields outside the particle boundaries. Such methods have been applied to simulations of particle-laden isotropic homogeneous flows^[5-6], pipe flows^[7], channel flows^[8-14], duct flows^[15-17], and Couette flows^[18-19].

The effects of finite-size neutrally buoyant particles on turbulent channel flows were investigated by many groups such as Shao et al.^[9], Picano et al.^[10], Wang et al.^[11], and Yu et al.^[12]. It was commonly observed that, in the near-wall region, the particles enhanced the transverse and spanwise root-mean-square (RMS) velocity fluctuations, but reduced the maximum streamwise RMS velocity. In the central region, the particle effects were opposite to those in the near-wall region. In addition, the flow drag was found to be enhanced for particle volume fractions of 1% and 10% for all simulations. In the simulations where the pressure gradient was fixed, the friction Reynolds number was normally chosen as 180, and the effects of the Reynolds number on the particle-turbulence interactions were not examined. The aim of the present study is to report our results on the effects of finite-size neutrally buoyant particles on the mean velocity, RMS velocities, the probability density function (PDF) of the velocity, and the vortex structures of the turbulent channel flow at the higher friction Reynolds number of 395. Further, the effects of the Reynolds number will be examined qualitatively, by comparing the results with those at Re_τ=180. We note that Re_τ=395 is a typical and widely used Reynolds number for the DNS of single-phase turbulent channel flows, similar to Re_τ=180, which is the reason why we choose it for our simulations.

2 Numerical model 2.1 Flow model

A schematic diagram for the channel flow studied is shown in Fig. 1. We let the x-axis represent the streamwise direction, the y-axis the wall-normal direction, and the z-axis the spanwise direction. The corresponding velocity components in the x-, y-, and z-directions are u, v, and w, respectively. The no-slip velocity boundary condition is imposed on the channel walls, and the periodic boundary condition is imposed in the streamwise and spanwise directions. We denote the half-width of the channel as H. In the present study, the computational domain is [0, 4H]×[-H, H] ×[0, 2H].

We take H as the characteristic length and the friction velocity as the characteristic velocity for the non-dimensionalization scheme. The friction velocity is defined as , where τ_w is the mean shear stress on the walls, and ρ_f is the fluid density. Thus, the Reynolds number is defined as Re_τ = u_τH/ν, where ν is the fluid kinematic viscosity. The pressure gradient is kept constant in our simulations, i.e., , and its dimensionless value normalized by ρ_fu_τ²/H is 1. The pressure gradient is imposed everywhere in the computational domain, including the exterior and interior of the particles.

2.2 Direct-forcing fictitious domain (DFFD) method

A DFFD method is used for the interface-resolved simulations of particle-laden turbulent flows^[20]. The fictitious domain (FD) method for the particulate flows was originally proposed by Glowinski et al.^[21]. The key idea of this method is that the interior of the particles is filled with the fluids, and the inner fictitious fluids are enforced to satisfy the rigid body motion constraint through a pseudo-body force, which is introduced as a distributed Lagrange multiplier in the FD formulation^[21]. In the following, we describe the DFFD method briefly, and the reader is referred to Yu and Shao^[20] for the details.

For simplicity, we consider only one particle in the following exposition. The particle density, volume and moment of inertia tensor, translational velocity, angular velocity, and position are denoted by ρ_s, V_p, J, U, ω_p, and X, respectively. Let P(t) represent the solid domain and Ω the entire domain, including the interior and exterior of the solid body. By introducing the following scales for the non-dimensionalization, H for the length, u_τ for the velocity, H/u_τ for the time, ρ_fu_τ² for the pressure, and ρ_fu_τ²/H for the pseudo-body force, the dimensionless FD formulations for the incompressible fluid can be written as follows:

(1)

(2)

(3)

(4)

(5)

In the above equations, u represents the fluid velocity, p is the fluid pressure with the exclusion of the constant pressure gradient of 1, λ is the pseudo-body force defined in the solid domain P(t), e_x represents the unit vector pointing to the streamwise direction, r is the position vector with respect to the mass center of the particle, ρ_r is the particle-fluid density ratio defined by ρ_r = ρ_s/ρ_f, Fr denotes the Froude number defined here by Fr = gH/u_τ², V_p^* is the dimensionless particle volume defined by V_p^* = V _p/H³, g represents the gravitational acceleration whose magnitude is denoted by g, and J^* represents the dimensionless moment of the inertia tensor defined by J^* = J/(ρ_sH⁵).

A fractional-step time scheme is used to decouple the system (1)-(5) into the following two sub-problems.

Consider the fluid sub-problem for u^* and p,

(6)

(7)

A finite-difference-based projection method on a homogeneous half-staggered grid is used for the solution of the above fluid sub-problem. All spatial derivatives are discretized with the second-order central difference scheme.

Consider the particle sub-problem for U^{n + 1}, ω_p^{n + 1}, λ^{n + 1}, and u^{n + 1},

(8)

(9)

Note that the above equations have been reformulated to ensure that all the right-hand side terms are known quantities. Consequently, the particle velocities U^{n + 1} and ω_p^{n + 1} are obtained without the iteration. Then, the pseudo-body force λ^{n + 1} defined at the Lagrangian nodes is determined by

(10)

Finally, the fluid velocity u^{n + 1} at the Eulerian nodes is corrected by

(11)

In the above manipulations, the tri-linear function is used to transfer the fluid velocity from the Eulerian nodes to the Lagrangian nodes, and the pseudo-body force from the Lagrangian nodes to the Eulerian nodes.

2.3 Collision model

A particle-particle collision model is required to prevent the mutual penetration of particles. We assume that a repulsive force is activated when the gap distance between two particles is smaller than a critical value. The force has the following form:

(12)

where F_ij, d_ij, and n_ij are the repulsive force acting on the particle j from the particle i, the gap distance, and the unit normal vector pointing from the center of the particle i to that of the particle j, respectively. Here, d_c represents a cut-off distance, the repulsive force is activated when d_ij < d_c, and F₀ is the magnitude of the force at contact. We set d_c=h, where h is the fluid mesh size, and F₀=10³. The motions of the particles due to the collision force (12) and due to the hydrodynamic force (8)-(11) are considered separately with a fractional step scheme. The time step for the collision model is set to be one-tenth of the latter (i.e., Δ t/10) to circumvent the stiffness problem arising from the explicit integration scheme with a large value of F₀, as suggested by Glowinski et al.^[21]. The collision between a particle and a wall is treated similarly with the coefficient F₀ in (12) doubled.

2.4 Validation

The accuracy of our DFFD code for the single-phase turbulence and good mesh-convergence of the turbulence statistics for the particle-laden flows (even for the coarse mesh of 6.4 points per one particle diameter) at Re_τ=180 were demonstrated by Yu et al.^[22]. Our results on the turbulent channel flow laden with finite-size neutrally buoyant particles were compared with lattice-Boltzmann simulations using interpolated bounce back at the fluid-solid interfaces, and the two completely different numerical approaches yielded quantitatively similar results in general^[11]. Here, we compare our results on the mean velocity and RMS velocities for the single-phase turbulent channel flow at Re_τ=395 with those of Iwamoto et al.^[23] from the pseudo-spectral simulations in Fig. 2, and good agreement can be seen.

3 Results and discussion

Throughout this study, the particle-fluid density ratio is ρ_r=1.0. We mainly report the results at Re_τ=395, but some results at Re_τ=180 are also presented for comparison. The particle radii normalized by the channel half-width are a/H = 0.05 and 0.025, corresponding to a⁺=au_τ/ν=19.75 and 9.875 (i.e., normalized by the wall unit), respectively. The particle volume fraction ϕ ranges from 0.098% to 7.07%. The parameter settings for the particle-laden cases at Re_τ=395 are presented in Table 1.

A homogeneous grid is used in our simulations. The grid number at Re_τ=395 is 512×256×256, corresponding to the mesh size h=H/128, and 12.8 points per particle diameter for a/H = 0.05. The grid number at Re_τ=180 is 256×128×128. It was shown that the coarse mesh of 6.4 points per one particle diameter could achieve good mesh-independent fluid statistics for the turbulent channel flow with a large particle-fluid density ratio^[22]. Thus, the results with the current mesh resolution should be reasonably accurate. The dimensionless time step is 0.000 1 at Re_τ=395 and 0.000 2 at Re_τ=180. The flow statistics are obtained from the averaging of the data in the real fluid domain over a period of typically 50 non-dimensional time units after the statistically steady state is reached.

Figure 3 shows the mean velocity profiles at Re_τ=395 for a/H=0.05 and 0.025. The relative values with respect to the particle-free case are also shown, to clarify the differences between the particle-laden cases at low volume fractions and the single-phase case. Because the pressure gradient is fixed in our simulations, a larger flow rate represents a lower flow friction. From Fig. 3(c), the flow rate in case of a/H=0.05 and ϕ=0.098% is slightly larger than that of the single-phase flow, indicating a slight drag reduction for this case. For all other cases, drag-enhancement is observed. Whether the spherical particles can produce drag-reduction has been debated for decades^[24-25]. The experimental data were contradictory^[24], and the numerical results of Zhao et al.^[25] indicated that the point particles could induce pronounced drag-reduction. Here, the drag-reduction by the finite-size particles is observed for the first time. Nevertheless, the amount of drag-reduction is negligibly small. Further studies are required to elucidate whether the finite-size spherical particles can produce pronounced drag-reduction. The flow drag for smaller particles at the same relatively large particle volume fraction is larger, as implied in Fig. 3 and also observed previously^[9]. The reason might be that the total surface area of smaller particles at the same particle volume fraction is larger, and then the viscous dissipation is larger, resulting in the larger flow friction.

As mentioned in the introduction, it was observed previously^[9-12] that the presence of particles decreased the maximum streamwise RMS velocity near the wall, while it increased the transverse and spanwise RMS velocities in the vicinity of the wall. In the center region, the particle effects were opposite to those in the near-wall region. The RMS velocity components for Re_τ=180 and a/H=0.05 are shown in Fig. 4, and qualitatively the same results can be clearly seen.

The RMS velocity components at Re_τ=395 are plotted in Fig. 5. The reduction in the maximum streamwise RMS velocity and the enhancement in the wall-normal and spanwise RMS velocities in the near-wall region (actually including the streamwise RMS velocity in the region very close to the wall) can be observed. Therefore, the particle effects on the RMS velocities are similar for different Reynolds numbers. However, the comparison between the RMS velocities at Re_τ=180 in Fig. 4 and those at Re_τ=395 in Fig. 5 indicates that the particle effects are weakened as the Reynolds number increases.

Typical vortex structures at Re_τ=395 for the particle-free and particle-laden cases are compared in Fig. 6. The vortex structures are identified with the λ_ci criterion^[26]. From Fig. 6, the particles weaken the large-scale vortices, which is expected to be responsible for the reduction in the maximum streamwise RMS velocity, and induce smaller-scale vortices, which causes the enhancements in all RMS velocity components near the wall. The weakening of large-scale vortices is probably a result of the particle-induced viscous dissipation and the disturbance from the particles. In general, the vortices of the single-phase turbulence at Re_τ=395 have smaller sizes and larger vorticity, and are more densely distributed in the channel, compared with the results at Re_τ=180. Therefore, the particles can affect more strongly the vortex structures and thereby the fluid statistics at Re_τ=180 than those at Re_τ=395.

Dinavahi et al.^[27] found that the PDFs of the fluctuating velocity components (beyond the buffer region) were roughly independent of the distance from the channel wall and the Reynolds number for the single-phase turbulent channel flow. To examine whether such universality holds for the particle-laden turbulent channel flow, the PDFs of fluctuating velocity components normalized by the RMS values at different y-positions are calculated and plotted in Fig. 7. The results for ϕ=2.36%, 7.07% are compared with those for the particle-free case and the Gaussian distribution. From Fig. 7, the normalized PDFs of the streamwise and wall-normal fluctuating velocities in and above the buffer region are different, but largely independent of the y-position above the buffer region. The normalized PDFs of the spanwise fluctuating velocity do not even depend on whether the y-position is in the buffer layer or not. Figure 7 further shows that the normalized PDFs are generally independent of the volume fraction and the particle size. Because the normalized PDFs for the particle-laden cases are close to those for the single-phase case, and the PDFs of the single-phase flow are independent of the Reynolds number^[27], it can be concluded that the normalized PDFs for the particle-laden cases are also independent of the Reynolds number.

4 Conclusions

We have investigated the effects of finite-size neutrally buoyant particles on turbulent channel flows for Re_τ=395 with a DFFD method. Our results show qualitatively the same particle effects at Re_τ=395 and 180, including the general drag-enhancement, the reduction in the maximum streamwise RMS velocity, and the enhancement in the wall-normal and spanwise RMS velocities in the vicinity of the wall. However, a very small drag-reduction is observed for a low particle volume at Re_τ=395, and our results indicate that particle effects on the mean and RMS velocities at Re_τ=395 are significantly smaller than those at Re_τ=180. The effects of the particles on the fluid velocity PDFs normalized with the RMS velocities are small, irrespective of the particle size, the particle volume fraction, and the Reynolds number.

The present study of the effects of the Reynolds number on the particle-turbulence interactions is preliminary. The effects of the Reynolds number on the particle statistics, and the quantitative effect of the Reynolds number, e.g., scalings, are good subjects for future studies.

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made.

The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.

To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.