Effect of channel dimensions and Reynolds numbers on the turbulence modulation for particle-laden turbulent channel flows

TL;DR

This study uses large eddy simulations to analyze how channel size and Reynolds number influence turbulence modulation in particle-laden turbulent flows, revealing key parameters affecting turbulence attenuation and structure.

Contribution

It identifies system size and fluid bulk Reynolds number as the primary factors influencing turbulence modulation in particle-laden flows, providing scaling insights for industrial design.

Findings

Turbulence attenuation varies with system size and Reynolds number.

Streamwise structures lengthen and become fewer with increased system size.

High Reynolds numbers produce smaller, thinner streaks that are more closely spaced.

Abstract

The addition of particles to turbulent flows changes the underlying mechanism of turbulence and leads to turbulence modulation. Different temporal and spatial scales for both phases make it challenging to understand turbulence modulation via one parameter. The important parameters are particle Stokes number, mass loading, particle Reynolds number, fluid bulk Reynolds number, etc., that act together and affect the fluid phase turbulence intensities. In the present study, we have carried out the large eddy simulations for different system sizes (2{\delta}/dp = 54, 81, and 117) and fluid bulk Reynolds numbers (Re_b = 5600 and 13750) to quantify the extent of turbulence attenuation. Here, {\delta} is the half-channel width, dp is the particle diameter, and Re_b is the fluid Reynolds number based on the fluid bulk velocity and channel width. The point particles are tracked with the Lagrangian approach. The scaling analysis of the feedback force shows that system size and fluid bulk Reynolds number are the two crucial parameters that affect the turbulence modulation more significantly than the other. The streamwise turbulent structures are observed to become lengthier and fewer with an increase in system size for the same volume fraction and fixed bulk Reynolds number. However, the streamwise high-speed streaks are smaller, thinner, and closely spaced for higher Reynolds numbers than the lower ones for the same volume fraction. In particle statistics, it is observed that the scaled particle fluctuations increase with the increase in system size while keeping the Reynolds number fixed. However, the scaled particle fluctuations decrease with the increase in fluid bulk Reynolds number for the same volume fraction and fixed system size. The present study highlights the scaling issue for designing industrial equipment for particle-laden turbulent flows.

Tables

Table 1. Table 1 : R ​ e b 𝑅 subscript 𝑒 𝑏 Re_{b} is the fluid bulk Reynolds number, δ 𝛿 \delta is half channel width, d p subscript 𝑑 𝑝 d_{p} is the particle diameter, S ​ t 𝑆 𝑡 St is the Stokes number based on the integral time scale of fluid ( 2 ​ δ / u ¯ 2 𝛿 ¯ 𝑢 2\delta/\bar{u} ), and S ​ t v 𝑆 subscript 𝑡 𝑣 St_{v} is the Stokes number based on the viscous time scale.

Case	$R{e}_b$	$2\delta /d_p$	$St$	$S{t}_v$
$A_1$	5600	54	179.53	4108.8
$B_1$	5600	81	74.45	1703.21
$C_1$	5600	117	38.71	885.27
$B_2$	13750	81	73.12	3304.4
$C_2$	13750	117	38.02	1720.8

Equations

$$\frac{\partial {\tilde{u}}_i}{\partial x_i}=0$$

$$\begin{gathered} \frac{\partial {\tilde{u}}_i}{\partial t}+\frac{\partial {\tilde{u}}_i{\tilde{u}}_j}{\partial x_j}=-\frac{1}{{\rho }_f}\frac{\partial \tilde{p}}{\partial x_i}+\nu \frac{{\partial }^2{\tilde{u}}_i}{\partial x_j\partial x_j} \\ +\frac{\partial ({\tilde{u}}_i{\tilde{u}}_j-\widetilde{u_i{u}_j})}{\partial x_j}+\frac{\widetilde{f_i}}{{\rho }_f}. \end{gathered}$$

$$\widetilde{f_i}=-\sum (F_{i,I}^d+F_{i,I}^l)$$

$${\nu }_t=C\tilde{\Delta }\sqrt{k_{SGS}}.$$

$$\frac{\partial k_{SGS}}{\partial t}+\frac{\partial (k_{SGS}{\tilde{u}}_i)}{\partial x_j}=\frac{\partial }{\partial x_j}\left[(\nu +{\nu }_t)\frac{\partial k_{SGS}}{\partial x_j}\right]-ϵ-{\tau }_{ij}\frac{\partial {\tilde{u}}_i}{\partial x_j}.$$

$${\tau }_{ij}-\frac{2}{3}{\delta }_{ij}{k}_{SGS}=-2C_1\tilde{\Delta }{k}_{SGS}^{1/2}{\tilde{S}}_{ij}.$$

$$T_{ij}-\frac{2}{3}{\delta }_{ij}{K}_{SGS}=-2C_2\tilde{\tilde{\Delta }}{K}_{SGS}^{1/2}{\tilde{\tilde{S}}}_{ij}.$$

$$L_{ij}=T_{ij}-\widetilde{{\tau }_{ij}}=\widetilde{\tilde{u_i}\tilde{u_j}}-\widetilde{\tilde{u_i}}\widetilde{\tilde{u_i}},$$

$$K_{SGS}=\widetilde{k_{SGS}}+\frac{1}{2}{L}_{ii}.$$

$$L_{ij}^d=L_{ij}-\frac{1}{3}{L}_{kk}{\delta }_{ij}={\alpha }_{ij}C-\widetilde{{\beta }_{ij}C},$$

$$M_{ij}={\alpha }_{ij}-\widetilde{{\beta }_{ij}}.$$

$$C=\frac{1}{2}\frac{⟨M_{ij}{L}_{ij}⟩}{⟨M_{kl}{M}_{kl}⟩}.$$

$$m_p\frac{d{v}_{i,I}}{dt}=F_{i,I}^D+F_{i,I}^L+\sum _{I\ne J}{F}_{i,IJ}+F_{i,Iw}+m_{p}g.$$

$$F_{i,I}^D=3\pi \mu d_p({\tilde{u}}_{i,I}(x,t)-v_{i,I})(1+0.15R{e}_p^{0.687})$$

$${\mathbf{F}}_{nij}=(-k_n{\delta }_n^{3/2}-{\eta }_{nj}G.n)\mathbf{n},$$

$${\mathbf{F}}_{tij}=-k_t{\delta }_t-{\eta }_{tj}{\mathbf{G}}_{ct}.$$

$$⟨\star {⟩}_s=\frac{1}{\delta }{\int }_0^{\delta }dy⟨\star ⟩.$$

$$-\frac{1}{{\rho }_f}\frac{\partial \tilde{P}}{\partial x}-\frac{\partial (\overline{\widetilde{u_x^{{}^{\prime }}}\widetilde{u_y^{{}^{\prime }}}})}{\partial y}+(\nu +{\nu }_t)\frac{{\partial }^2({\tilde{U}}_x)}{\partial y\partial y}-\overline{\frac{{\rho }_{p}f{\varphi }_l}{{\rho }_f{\tau }_p}({\tilde{u}}_x-{\tilde{v}}_x)}=0,$$

$$-\frac{1}{{\rho }_f}\frac{\partial \tilde{P}}{\partial y}-\frac{\partial (\overline{\widetilde{u_y^{{}^{\prime }}}\widetilde{u_y^{{}^{\prime }}}})}{\partial y}-\overline{\frac{{\rho }_{p}f{\varphi }_l}{{\rho }_f{\tau }_p}({\tilde{u}}_y-{\tilde{v}}_y)}=0.$$

$$\begin{gathered} (\nu +{\nu }_t)\frac{\partial }{\partial y}({\tilde{U}}_x\frac{\partial \widetilde{U_x}}{\partial y})-\frac{\partial (\overline{\widetilde{u_x^{{}^{\prime }}}\widetilde{u_y^{{}^{\prime }}}}{\tilde{U}}_x)}{\partial y}+\overline{\widetilde{u_x^{{}^{\prime }}}\widetilde{u_y^{{}^{\prime }}}}\frac{\partial ({\tilde{U}}_x)}{\partial y} \\ -(\nu +{\nu }_t)\frac{\partial {\tilde{U}}_x}{\partial y}\frac{\partial {\tilde{U}}_x}{\partial y}-{\tilde{U}}_x\overline{\frac{{\rho }_{p}f\varphi }{{\rho }_f{\tau }_p}({\tilde{u}}_x-{\tilde{v}}_x)}-{\tilde{U}}_x\frac{1}{\rho }\frac{\partial \tilde{P}}{\partial x}=0. \end{gathered}$$

$$\begin{gathered} ⟨\overline{\widetilde{u_x^{{}^{\prime }}}\widetilde{u_y^{{}^{\prime }}}}\frac{\partial ({\tilde{U}}_x)}{\partial y}{⟩}_s-(\nu +{\nu }_t)⟨\frac{\partial {\tilde{U}}_x}{\partial y}\frac{\partial {\tilde{U}}_x}{\partial y}{⟩}_s \\ -⟨{\tilde{U}}_x\overline{\frac{{\rho }_{p}f{\varphi }_c}{{\rho }_f{\tau }_p}({\tilde{u}}_x-{\tilde{v}}_x)}{⟩}_s-{\tilde{U}}_x\frac{1}{{\rho }_f}\frac{\partial \tilde{P}}{\partial x}=0 \end{gathered}$$

$$F_p=\frac{{\rho }_{p}f{\varphi }_l}{{\rho }_f{\tau }_p}({\tilde{u}}_x-{\tilde{v}}_x)\frac{2\delta }{{\overline{u}}^2},$$

$$F_p=\frac{{\rho }_{p}f{\varphi }_l}{{\rho }_f({\rho }_p{d}_p^2/18\mu )}({\tilde{u}}_x-{\tilde{v}}_x)\frac{2\delta }{{\overline{u}}^2},$$

$$F_p=18f{\varphi }_l\frac{R{e}_p}{R{e}_b^2}{\left(\frac{2\delta }{d_p}\right)}^3.$$

$$F_p\propto {\delta }^3.$$

$$F_p\propto \frac{1}{R{e}_b^2}.$$

$$r_{ij}(r,t)=\overline{u_i^{\prime }(x,t)u_j^{\prime }(x+r,t)}.$$

$$R_{ij}(r,t)=\frac{\overline{u_i^{\prime }(x,t)u_j^{\prime }(x+r,t)}}{\overline{u_i^{\prime }(x,t)u_j^{\prime }(x,t)}}.$$

$$D_{ij}=\frac{⟨u_i^{\prime }(0)u_j^{\prime }(0)⟩}{{\tau }_p^2}\underset{0}{\overset{\infty }{\int }}d{t}^{\prime }{R}_{ij}.$$

$$D_{ij}=\frac{⟨u_i^{\prime }(0)u_j^{\prime }(0)⟩}{{\tau }_p^2}{\tau }_I.$$

Full text

particle-laden flows, LES, turbulence modulation,

††preprint: AIP/123-QED

I Introduction

The Particle laden turbulent flows are ubiquitous in different industrial processes such as combustion, coating, fluidization, pneumatic transport of solids, electrostatic precipitators, and natural processes like raindrop formation, dust storms, etc. The complex interactions between the gas and solid phases lead to the preferential accumulation of high inertia particles at low vorticity zones, crossing trajectories, turbophoresis induced inhomogeneous distribution of particles, and fluid phase turbulence modulation. Turbulence modulation is an important area of investigation as the extent of turbulence modulation directly impacts the heat and mass transfer properties. Turbulence modulation depends on many parameters such as Stokes number, particle Reynolds number, fluid Reynolds number, the ratio of integral length scale to the particle diameter, solid mass loading, etc. In an early work, Gore and Crowe (1989) divided the augmentation and attenuation region based on the integral length scale ratio to the particle diameter ($L/d_{p}$). The authors found that if the ratio is less than 0.1, the attenuation will be observed; otherwise, augmentation will be observed. Tanaka and Eaton (2008) presented two dimensionless numbers based on the Stokes number and particle Reynolds number ($Re_{p}$), and differentiated the regions of augmentation and attenuation based on these two dimensionless numbers. They found that region of attenuation exist between two augmentation regions. Yu et al. (2021a) have proposed the criteria to predict the attenuation or augmentation for turbulent kinetic energy. The authors have performed fully resolved simulations for an upward channel for two bulk Reynolds numbers ($Re_{b}$) of 5746 and 12000, $\rho_{p}/\rho_{f}=2-100$, and $\phi=0.003-0.0236$. Here, $\phi$ is the solid volume fraction, $\rho_{p}$ is the particle density, and $\rho_{f}$ is the fluid density. For $Re_{p}<50$, $\rho_{p}/\rho_{f}=2$, and $2\delta/d_{p}=20$, an increase in bulk Reynolds number from 5746 to 12000 leads to less attenuation for $\phi=0.0236$ (Fig. 11 of their paper). The authors observed that at low $Re_{p}$, the turbulence attenuation is observed across the channel; at intermediate $Re_{p}$, the turbulence intensity is enhanced in the channel center and diminished in the near wall. However, at sufficiently high $Re_{p}$, the turbulence intensity is increased across the channel width. Yu et al. (2017) have simulated the fully resolved DNS at friction Reynolds number of 180, $\phi=0.84\%$, and $\rho_{p}/\rho_{f}=1-104.2$. The effect of density ratio on the turbulence modulation in a turbulent channel flow with constant pressure gradient and without considering the gravity effect is studied. The drag modification with varying density ratio is not monotonic, and authors found large drag at a density ratio of 10.42 compared to unity and 104.2. It is mentioned that the reduction in Reynolds stress results in drag reduction, but it is counteracted by the increase in particle-induced stress and leads to an increase in overall dragYu et al. (2017); Costa et al. (2021).

In their experimental work, Kulick et al. (1994) found that turbulence attenuation increases with increased Stokes number, mass loading, and distance from the wall. The authors studied the turbulence modulation at $Re_{\tau}=644$ based on half-channel width and friction velocity. They didn’t observe the attenuation in fluid fluctuations for low Stokes numbers (fluid time scale based on $k/\epsilon$ at the centerline) of 0.57 (glass particles). However, a significant decrease in the streamwise fluid fluctuations was observed for Stokes number of 3 (copper particles) at a mass loading of 0.8. Li et al. (2001) carried out DNS studies at low $Re_{\tau}=125$ for the vertical channel flows and studied the effect of mass loading, Stokes number, and density ratio. A complete turbulence collapse was observed for the high mass loading, $\phi_{m}=2$. Yamamoto et al. (2001) performed large eddy simulations (LES) for a vertical channel at Reynolds number ($Re_{\tau}$) of 644 and emphasized the importance of particle-particle collisions even for low volume fraction cases. The authors observed very small attenuation even at mass loading of one for the copper particle of $St=\tau_{p}U_{cl}/\delta$ = 70 where $\tau_{p}$ is the particle relaxation time, $\delta$ is the half-channel width, and $U_{cl}$ is the channel centerline velocity. Duque-Daza et al. (2021) performed the LES of turbulent channel flows for three different Reynolds numbers of $Re_{\tau}=180$, 365, and 950. For mass loading of 2.96, a complete turbulence collapse was observed at $Re_{\tau}=180$, while a negligible turbulence modification was observed at $Re_{\tau}=950$. The Stokes number can be expressed as $St=\tau_{p}/\tau_{k}=(1/18)(d_{p}/\eta)^{2}(\rho_{p}/\rho_{f})$. Here, $\tau_{k}$ is the Kolmogorov velocity scale, $\eta$ is the Kolmogorov length scale, and $\tau_{p}$ is the particle relaxation time. Thus, the Stokes number can be modified by varying the density ratio or changing particle diameter. The effect of both the parameters was investigated by Shen et al. (2022) for homogenous isotropic turbulence with fully resolved simulations. In their study, it was observed that particles always lead to turbulence attenuation; attenuation was found to be larger for high-density particles when the particle diameters were the same. However, the attenuation is smaller for the larger particles for the fixed density ratio due to vortex shedding behind the particles when the particle Reynolds number is significantly larger. Lee and Lee (2015) performed DNS to examine the effect of Stokes number on the turbulence modulation for the channel flow. The authors observed that the smaller particles ( $St^{+}=0.5$) act as the energy source for the high-low speed streaks that may increase the instability and gives rise to the new quasi-streamwise vortices. However, the turbulence was decreased for the high $St^{+}$ particles. An increase in the number of quasi-streamwise vortices was observed at $St^{+}=0.5$, while the number is less for the high Stokes numbers.

The effect of particle loading on the turbulent structure, viz. streak spacing in turbulent channel flows, has been described by different authors Zhao et al. (2010); Zhou et al. (2020); Duque-Daza et al. (2021). Zhou et al. (2020) observed a non-monotonic behavior of streamwise fluid fluctuations with an increase in mass loading. However, the wall-normal and spanwise fluid fluctuations decreased with an increase in particle mass loading. It was observed that the streamwise fluctuations increased up to a mass loading of 0.75, and a further increase in mass loading to 0.96 caused a decrease in streamwise fluid fluctuations. The authors presented two phenomena related to the number and spacing of the near-wall streaks and modifications of the strength of vortices. First, the near-wall vortices become fewer and wider, and the spacing between the streaks increases with an increase in mass loading. This results in a decrease in the streamwise fluid fluctuations. Second, the vortices and streaks aligned in the streamwise direction, and enhanced the strength of the streaks, increasing the streamwise fluctuations beyond the mass loading of 0.75. Zhao et al. (2010) and White and Mungal (2008) mentioned that the particles affect the regeneration cycle of near-wall turbulence by reducing the strength and the number of quasi-streamwise vortices. This leads to wider and more stable near-wall streaks. Dritselis and Vlachos (2008) concluded that the particles provide torque in the opposite direction to that of the mean vortex, which leads to the weaker and wider mean fluid structure. Ghosh and Goswami (2022) discussed the role of inter-particle collisions in breaking up the high-low speed streaks due to particles in the context of particle-laden shear flows. The turbulence modulation and particle clustering phenomena have been reviewed by Balachandar and Eaton (2010); Brandt and Coletti (2022). A number of studies on particle-laden turbulent channel flows have also been summarized by Muramulla et al. (2020), Wang et al. (2020), and Rohilla et al. (2022).

Tanaka and Eaton (2008) have proposed two dimensionless numbers (particle momentum numbers) to classify the regimes of turbulence modulation, which are based on the particle Reynolds number and Stokes number. One of the particle momentum numbers based on the Stokes number is given as, $Pa_{St}=\frac{1}{54\sqrt{2}}\frac{Re_{b}^{2}}{St_{K}^{1/2}}\left(\frac{\rho_{p}}{\rho_{f}}\right)^{3/2}\left(\frac{d_{p}}{L}\right)^{3}$, when $v_{k}\sim v_{rel}$ in the limit of $Re_{p}\sim d_{p}/\eta$. Here, $St_{K}$ is the Stokes number based on the Kolmogorov time scale, $v_{k}$ and $\eta$ are the Kolmogorov velocity and length scales, and $v_{rel}$ represents the relative velocity between the fluid and particle phases. Another dimensionless number based on the particle Reynolds number is defined as $Pa_{Re}=\frac{1}{18}\frac{\rho_{p}}{\rho_{f}}\left(\frac{d_{p}}{L}\right)^{3}\frac{Re_{b}^{2}}{Re_{p}}$. The particle Reynolds number is not an independent variable. However, other parameters such as density ratio, bulk Reynolds number, and the ratio of particle size to integral length scale can be input parameters. The effect of mass loading or the density ratio on the turbulence modulation for the channel flows Kulick et al. (1994); Li et al. (2001); Wu et al. (2006); Lee and Lee (2015); Capecelatro, Desjardins, and Fox (2018); Fong et al. (2019); Muramulla et al. (2020); Zhou et al. (2020), homogenous isotropic turbulence Ten Cate et al. (2004); Yeo et al. (2010); Shen et al. (2022), homogenous shear Ahmed and Elghobashi (2000); Yousefi, Ardekani, and Brandt (2020), decaying turbulence Lucci, Ferrante, and Elghobashi (2010); Gao, Li, and Wang (2013); Schneiders, Meinke, and Schröder (2017); Luo et al. (2017), Couette flow Richter and Sullivan (2013, 2014); Richter (2015); Ghosh and Goswami (2022) have been reported. The different particle-laden experimental works are summarized by Wu et al. (2006) and Fong et al. (2019). Effect of channel dimensions (ratio of the integral length scale to the particle diameter) and Reynolds number on the turbulence modulation in the limit of high Stokes number and low solid volume loading is scarce in the literature Yu et al. (2017, 2021b). In the case of LES, there are limited studies where turbulence modulation has been explored with the variation of volume fraction/mass loading and Reynolds number Yamamoto et al. (2001); Duque-Daza et al. (2021); Rohilla et al. (2022). In this work, we have studied the effect of $L/d_{p}$ and Reynolds numbers on turbulence modulation and found that these are the two important parameters that quantify the attenuation level relative to volume fraction. The present study highlights the important scaling issue which should be dealt while comparing the simulation and experimental results.

Paper outline: In the second section, the governing equations for fluid and particle phases are discussed briefly, along with the simulation parameters used in the present study. The third section presents the turbulence modulation of fluid fluctuations and particle properties due to the variation in volume fractions, system sizes, and Reynolds numbers. The conclusions are discussed in the last section.

II Governing equations

II.1 Fluid Phase Equations

The incompressible fluid phase is described in the Eulerian framework using Navier-Stokes equations. In the LES approach, the instantaneous mass and momentum equations are defined using the filtered quantities written as,

[TABLE]

and,

[TABLE]

Where $\widetilde{u}_{i}$ and $\widetilde{p}$ are the filtered velocity and pressure, respectively, $\nu$ is the kinematic viscosity, and $\rho_{f}$ is the fluid density. The ($\tilde{.}$) denotes the filtered quantity throughout the manuscript. $\widetilde{f_{i}}$ is the feedback force per unit volume caused by the dispersed phase. The feedback force includes the drag and the lift on the particles, and is written as,

[TABLE]

Here, $F_{i,I}^{d}$ and $F_{i,I}^{l}$ are the drag and lift forces on the $I^{th}$ particle. The third term on the right-hand side of Eqn. (2) is referred to as the subgrid-scale (SGS) stress term, which requires closure. The dynamic one-equation model has been used as a fluid SGS model. A correction on the lift force term is applied following Mei (1992) and Loth and Dorgan (2009). The importance of the inclusion of lift force is discussed by Costa et al. (2020). In the present methodology, a finite volume based opensource software CFDEM Kloss et al. (2012); El-Emam et al. (2021) has been used. A second-order backward scheme for time derivative and a second-order Gauss linear differencing scheme for pressure gradient, viscous diffusion, and advection terms have been adopted.

II.1.1 Dynamic one-equation model

In the dynamic one-eqaution model, the subgrid scale stress ($\tau_{ij}$) is expressed in terms of the eddy viscosity ($\nu_{t}$) and filtered strain rate ($\widetilde{S}_{ij}$) as $\tau_{ij}=-2\nu_{t}\widetilde{S}_{ij}$. The eddy viscosity is calculated as,

[TABLE]

Here, $k_{SGS}=\frac{1}{2}[\widetilde{u_{i}u_{i}}-\widetilde{u}_{i}\widetilde{u}_{i}]$ is the subgrid scale kinetic energy, $\widetilde{\Delta}$ is the cube root of grid volume, and $C$ is a parameter calculated dynamically during the simulation. The $k_{SGS}$ is calculated using the following equation (Eqn. 5).

[TABLE]

Where $\epsilon=C_{\epsilon}k_{SGS}^{3/2}/\widetilde{\Delta}$, and $C_{\epsilon}=1.048$ is the model coefficient. In this model, the model constant appearing in eddy viscosity ($C$ in Eqn. 4) is calculated dynamically Pomraning and Rutland (2002) which is explained as,

[TABLE]

Using the second or test filter, the subgrid stress at the test filter level can be written as,

[TABLE]

Where $K_{SGS}=\frac{1}{2}[\widetilde{\widetilde{u_{i}u_{i}}}-\widetilde{\widetilde{u}}_{i}\widetilde{\widetilde{u}}_{i}]$. The calculation of $C$ involves Germano’s identity Germano et al. (1991),

[TABLE]

and,

[TABLE]

The formulation involves two assumptions. First, the subgrid stress scales with filter size, then $C_{1}$ and $C_{2}$ are replaced by $C$. And second, $C$ is taken outside the integral in Eqn. 10 as the variation of $C$ is low in space. Using Eqns. 6 - 8, the deviatoric part of the $L_{ij}$ is written as,

[TABLE]

and,

[TABLE]

Where $\alpha_{ij}=-2\widetilde{\widetilde{\Delta}}K_{SGS}^{1/2}\widetilde{\widetilde{S_{ij}}}$ and $\beta_{ij}=-2\widetilde{\Delta}k_{SGS}^{1/2}\widetilde{S_{ij}}$. Then,

[TABLE]

Here, $C_{s}(=\sqrt{C})$ is the Smagorinsky coefficient. The angular brackets in Eqn. 12 denote the plane averaging. The present work uses a box filter as a test filter. This model does not require the Van-driest damping as in the Smagorinsky model to satisfy the wall scaling.

II.2 Particle phase description

The solid individual particle is tracked in the Lagrangian framework. The particle-particle and particle-wall collisions are considered to be elastic in nature. The particle motion of the dispersed phase is described using Newton’s second law of motion, which is defined as,

[TABLE]

Where $m_{p}$ is the mass of the particle, $v_{i,I}$ is the velocity of $I^{th}$ particle, $F_{i,I}^{D}$ and $F_{i,I}^{L}$ are the drag and lift force acting on the particle, respectively. The g is the gravitational acceleration, $F_{i,IJ}$ is the interaction force between $I^{th}$ and $J^{th}$ particle, and $F_{i,Iw}$ is the interaction force between $I^{th}$ particle and the wall. The drag force is calculated using the Schiller-Naumann correlation Naumann and Schiller (1935) given by Eqn. 14.

[TABLE]

A third-order accurate scheme has been used to interpolate the fluid velocity at the particle location to calculate the drag and the lift forces. The particle size is smaller than the grid size. A point particle approximation has been used for the present study Boivin, Simonin, and Squires (2000); Balachandar (2009); Balachandar and Eaton (2010); Zhao et al. (2010); Dritselis and Vlachos (2011); Mallouppas and van Wachem (2013); Marchioli (2017).

The soft-sphere (spring dashpot) Cundall and Strack (1979); Crowe et al. (1998) interaction model has been used to capture the particle-particle and particle-wall collisions. The particle collisions are important to capture even at a low solid volume fraction($\phi\sim 10^{-4}$) Yamamoto et al. (2001) . The collisions lead to momentum transfer to all the directions and affect the particle phase statistics. During the collision, the repulsive force on the particle depends on the extent of overlapping of the colliding particles. This is modeled using a spring-dashpot system where the dashpot includes the energy loss associated with the collision. The motion of particles is divided into normal and tangential components, along with overlapping distances and forces. Considering two colliding particles as $i$ and $j$, the normal ($\mathbf{F}_{nij}$) and tangential ($\mathbf{F}_{tij}$) components of the forces are modeled as,

[TABLE]

[TABLE]

Here, the subscripts ‘n’ and ‘t’ refer to the normal and tangential directions. $\delta$ is the overlap distance, $k$, and $\eta$ are the stiffness and damping coefficient, $\mathbf{G}$ is the relative velocity vector of $i^{th}$ particle relative to $j^{th}$ particle ($\mathbf{G=v_{i}-v_{j}}$), and $\mathbf{n}$ is the unit vector in the direction of the line connecting the centers of both the particles. $\mathbf{G}_{ct}$ is the slip velocity at the contact point, and defined as $\mathbf{G}_{ct}=\mathbf{G-(G.n)n}$.

II.3 Simulation parameters

Simulations are performed at two bulk Reynolds numbers ($Re_{b}=\bar{u}\times 2\delta/\nu$) and three different channel dimensions ($2\delta/d_{p}$). Where $\delta$ is the half-channel width, $\bar{u}$ is the average fluid velocity, and $d_{p}$ is the particle diameter. The bulk Reynolds numbers are 5600 and 13750, and three different channel dimensions used are 54, 81, and 117. Lengths of the vertical channel considered are $8\pi\delta$, $2\delta$, and $(4/3)\pi\delta$ in streamwise (x), wall-normal (y), and spanwise (z) directions, respectively. The finite volume method is applied to solve the filtered Navier-Stokes Eqns. 1 and 2. The fluid phase is resolved with $128*65*64$ and $192*99*96$ grid points ($x*y*z$) for bulk Reynolds numbers of 5600 and 13750, respectively. The grid resolution is $\triangle x^{+}=35$, $\triangle z^{+}=12$ for $Re_{b}=5600$, and $\triangle x^{+}=51$, $\triangle z^{+}=17$ for $Re_{b}=13750$ for the LES models. Earlier studiesWang and Squires (1996); Salmanzadeh et al. (2010) have also used similar grid resolution . Here, (+) symbol indicates the quantities normalized with viscous scales. In the wall-normal direction, the first grid point is such that $y^{+}$ is less than one. No-slip boundary condition is applied at the walls. The time step is taken such that the CFL number is nearly 0.2 for all the constant bulk flow rate simulations. The Dynamic one-equation (Eqns. 4 - 12) LES model is used to close the SGS stress term in the filtered Navier-stokes equation.

The particle diameter is considered to be $39\mu m$. The particle densities are considered as $2000kg/m^{3}$ and $800kg/m^{3}$. As the particle to fluid density ratio is higher, the Basset history and buoyancy effects are neglected in the particle’s equation of motion. Simulation parameters like Reynolds numbers, channel dimensions, and the particle Stokes numbers are provided in Table 1. The Stokes number based on the fluid integral length scale and viscous scales are denoted by $St=\tau_{p}/\tau_{f}$ and $St_{v}=\tau_{p}/\tau_{fv}$. Here, the particle relaxation time ($\tau_{p}$) is defined as, $\tau_{p}=\rho_{p}d_{p}^{2}/18\mu_{f}$, the fluid integral time scale is defined as, $\tau_{f}=2\delta/\overline{u}$, and viscous time scale is defined as $\tau_{fv}=\nu/u_{\tau}^{2}$ where $u_{\tau}$ is the friction velocity. Here, unladen friction velocity is taken to calculate $St_{v}$ in Table 1. Various authors (Yamamoto et al., 2001; Kuerten and Vreman, 2005; Kuerten, 2006; Marchioli, Salvetti, and Soldati, 2008; Dritselis and Vlachos, 2011; Zhou et al., 2020) have used time scale based on the friction velocity to define the Stokes number. However, it is to be noted that the turbulence modulation with the increase in particle mass loading leads to a modification in friction velocity, shown in Fig. 1(d) (a). For that reason, we have used $2\delta/\bar{u}$ as the fluid time scale.

III Results

III.1 Effect of system size

In this section, the simulations are performed to analyze the dynamics of fluid and solid phases with different channel widths for a range of solid volume fractions at $Re_{b}=5600$. The channel dimensions are considered such that the ratios of channel width ($2\delta$) and particle diameter ($d_{p}$) are 54, 81, and 117. Here, channel width ($2\delta$) varies for different channel configurations. However, the computational domain is kept as $8\pi\delta\times 2\delta\times(4/3)\pi\delta$ for all the cases. The comparison of current unladen LES with available literature Moser et al. (1999); Salmanzadeh et al. (2010); Li et al. (2012); Muramulla et al. (2020) is shown in the appendix for $Re_{b}=5600$ and 13750. In Fig. 1(d) (a), the fluid streamwise mean velocity is plotted for all three channels and two volume fractions, $\phi=2\times 10^{-4}$ and $5\times 10^{-4}$. It is observed that there is no variation in the mean profiles for the case $A_{1}$ when the volume fraction is increased from $\phi=2\times 10^{-4}$ to $5\times 10^{-4}$. For case $B_{1}$, there is a marginal decrease of about 5% in mean velocity in the buffer region as the volume fraction is increased from $\phi=2\times 10^{-4}$ to $5\times 10^{-4}$. However, the reduction in the mean velocity is more significant for the case $C_{1}$ at $\phi=5\times 10^{-4}$. The normalized fluid fluctuations as a function of the wall-normal distance are plotted in Fig. 1(d) (b-d). In the case of streamwise fluid fluctuations, there is only 4% and 9% decrease in the peak value as the volume fraction is increased from $\phi=2\times 10^{-4}$ to $5\times 10^{-4}$ for the case $A_{1}$ and $B_{1}$, respectively. However, the decrease is almost one order of magnitude for the case $C_{1}$ for the same increase in volume fraction (Fig. 1(d) (b)). In the case of cross-stream fluid velocity fluctuations, there is 13% and 34% decrease as the volume fraction is changed from $\phi=2\times 10^{-4}$ to $5\times 10^{-4}$ for the case $A_{1}$ and $B_{1}$, respectively. However, the decrease is of almost one order of magnitude for similar volume fraction change in case of the larger channel of $2\delta/d_{p}=117$ (case $C_{1}$), Fig. 1(d) (c). A similar behavior is observed for wall-normal fluid fluctuations, shown in Fig. 1(d) (d).

To quantify the variation of turbulence intensity, we compute channel averaged fluid fluctuations as defined by

[TABLE]

In Fig. 2(d) (a-d), the average fluid fluctuations across the channel width are plotted for a range of volume fractions. The average fluctuations predicted by different channels match at lower volume fractions. The fluid fluctuations decrease with an increase in volume fraction, and a complete turbulence suppression is observed at a higher particle volume loading for all three cases. The particle loading at which complete turbulence collapse occurs is called critical particle volume loading (CPVL). The attenuation observed in the fluid fluctuations is higher for the case $C_{1}$. The CPVL is $14\times 10^{-4}$, $7\times 10^{-4}$, and $5\times 10^{-4}$ for the cases $A_{1}$, $B_{1}$ and $C_{1}$, respectively. It is interesting to note that the fluid fluctuations decrease with an increase in system size for the same volume fraction and fixed Reynolds number, although, Stokes number is decreased (Table 1). Such an observation is in contrary to the earlier findingsLi et al. (2001); Dritselis and Vlachos (2011); Lee and Lee (2015); Muramulla et al. (2020); Kumaran et al. (2020) that the fluid fluctuations decrease with an increase in Stokes number for $St_{K}>1$. Here, $St_{K}$ is the Stokes number based on the Kolmogorov time scale. It is worth noting that, in the earlier studies Stokes number was changed by changing particle relaxation time for a fixed fluid time scale. However, in the present study, particle relaxation time is fixed and the fluid time scale changes with the system size. The results show that there is a reduction in the fluid fluctuations with the decrease in considered Stokes numbers. This raises the question whether the Stokes number is a correct parameter to represent the effect of particle inertia on the turbulence modulation for a range of system size and in the limit $St_{K}>1$. Addressing such a question is important if one is interested to quantify turbulence modulation in an unified manner as a function of fluid and particle inertia. Another interesting implication is that the understanding gained will help to compare the simulation results with the experiments.

III.2 Effect of Reynolds number of fluid phase

This section examines the turbulence modulation for two channel dimensions of $2\delta/d_{p}=81$ and 117. The unladen fluid statistics are validated against the DNS data of Moser et al. (1999), experimental data of Li et al. (2012), and LES data of Salmanzadeh et al. (2010) in the appendix, shown in Fig. 20(d). In Fig. 3(d), the fluid statistics for two volume fractions of $\phi=5\times 10^{-4}$ and $10^{-3}$ are presented for cases $B_{1}$ and $B_{2}$. In Fig. 3(d) (a), it is observed that there is no change in the streamwise mean velocity profiles as the volume fraction is increased from $\phi=5\times 10^{-4}$ to $10^{-3}$ for the case $B_{2}$. However, in case of $B_{1}$, there is a decrease of fluid streamwise mean velocity in the buffer region and an increase in the mean velocity near the channel center for a similar change in volume fraction (shown in Fig. 3(d) a). Fig. 3(d) shows that there is about a 5% decrease in the peak value of streamwise fluid fluctuations as the volume fraction is increased from $\phi=5\times 10^{-4}$ to $10^{-3}$ for the case $B_{2}$. However, a decrease of almost one order of magnitude in the peak value is observed for the case $B_{1}$ for the same variation of volume fraction. There is a decrease of nearly 14% and 11% in the peak value of the Reynolds stress and wall-normal fluctuations, respectively, as the volume fraction is increased from $\phi=5\times 10^{-4}$ to $10^{-3}$ for the $Re_{b}=13750$. However, in the case of $Re_{b}=5600$, there is a decrease of almost one order of magnitude in the peak value of the Reynolds stress and the wall-normal fluctuations.

Since the extent of variation of second-order velocity fluctuations is a function of wall-normal distance, we compute a channel averaged value of the moments to quantify the variation of the moments as a function of solid volume fraction. The average fluid fluctuations are plotted in Fig. 4(d) for both the Reynolds numbers of 5600 and 13750, and for channel dimensions of $2\delta/d_{p}=81$ and 117 over a range of particle volume fractions. The decrease in the fluid fluctuations is observed with an increase in particle volume loading. It is observed that there is almost a 40% decrease in the streamwise fluid fluctuations and almost 60% decrease in the cross-stream, wall-normal, and spanwise fluid fluctuations compared to unladen cases before the CPVL for both the channels and Reynolds numbers. However, for the constant $2\delta/d_{p}$, the CPVL for $Re_{b}=13750$ is significantly higher than that for the $Re_{b}=5600$.

To understand the dependence of turbulence modulation and its collapse on the channel dimension and fluid phase Reynolds number, we conduct momentum and energy balance for a range of solid volume fractions and discuss them in the following section.

III.3 Momentum and energy balance

This section discusses the terms of the momentum balance and energy budget for the mean fluid flow. The filtered momentum equations for streamwise and wall-normal components are written as,

[TABLE]

and

[TABLE]

Here, $\widetilde{p}$ and $\widetilde{u}_{i}$ are the fluid’s instantaneous filtered pressure and velocity, respectively. $\widetilde{P}$ and $\widetilde{U}_{x}$ are the mean filtered pressure and velocity, respectively. The filtered fluctuating velocity is denoted by $\widetilde{u^{\prime}}_{i}$. The last term in both the above equations shows the feedback force exerted by the particles, where $\rho_{p}$ and $\rho_{f}$ are the particle and fluid density, $f(=1+0.15Re_{p}^{0.687})$ is the drag factor, $\tau_{p}$ is the particle relaxation time, and $\phi_{l}$ is the local volume fraction of the particle. The filtered Reynolds stress is represented by $\rho_{f}\overline{\widetilde{u^{{}^{\prime}}_{i}}\widetilde{u^{{}^{\prime}}_{j}}}$. The terms appearing in the streamwise filtered momentum equation are plotted in Fig. 5(c) as a function of wall-normal positions. All the terms are scaled with $\bar{u}^{2}/2\delta$. The wall-normal distance ($y$) is normalized with channel width ($2\delta=h$). The terms of the momentum equation are plotted for different channel dimensions and volume fractions which present the effect of system size. Here, the derivative of Reynolds stress is denoted by $R_{12,y}$, diffusion term by $D_{m}$, pressure gradient by $\Pi_{m}$, and feedback force by $F_{p}$. For unladen cases, the terms for channels with different dimensions match with each other as expected (Fig. 5(c) a). Here, magnitude of $R_{12,y}$ is comparable to the $D_{m}$. However, the deviation is observed at $\phi=2\times 10^{-4}$ in Fig. 5(c) (b), and the peak values of $R_{12,y}$ and $D_{m}$ are smaller by approximately 20% for case $C_{1}$ compared to the case $A_{1}$ (Table 1). A further increase in volume fraction to $\phi=5\times 10^{-4}$ results a decrease in $D_{m}$, and $R_{12,y}$ almost becomes zero for case $C_{1}$ (higher $2\delta/d_{p}$) as shown in Fig. 5(c) (c). Fig. 5(c) (c) shows that $\langle\Pi_{m}\rangle_{s}$ decreases nearly 30%, while $\langle F_{p}\rangle_{s}$ increases almost by 400% for case $C_{1}$ compared to case $A_{1}$. Here, $\langle.\rangle_{s}$ denotes the averaging of a quantity across the channel width as defined by Eqn. 17. Here, a significant increase in the feedback force is observed with an increase in $2\delta/d_{p}$.

To quantify the effect of the Reynolds number, we have chosen a channel with $2\delta/d_{p}=81$, and the terms of the momentum equation are plotted for $Re_{b}=5600$ and 13750. In the case of the unladen flow, the peak location of $R_{12,y}$ and $D_{m}$ are closer to the wall, and magnitude of peak values are almost twice for $Re_{b}=13750$ as shown in Fig. 6(c) (a). It is interesting to note that the scaled pressure gradient is lowered by approximately 25% for the case $B_{2}$ ($Re_{b}=13750$) than case $B_{1}$ ($Re_{b}=5600$), Fig. 6(c) (a). A decrease of nearly 25% in $\langle\Pi_{m}\rangle_{s}$ and $\langle D_{m}\rangle_{s}$ is observed, while there is no variation of $\langle R_{12,y}\rangle_{s}$ in case of $B_{2}$ compared to $B_{1}$. The terms are shown for $\phi=5\times 10^{-4}$ and $10^{-3}$ in Fig. 6(c) (b and c). $R_{12,y}$ and $D_{m}$ decrease with an increase in volume fraction. However, the extent of the decrease is much lower for $Re_{b}=$13750 than $Re_{b}=5600$. $R_{12,y}$ is almost zero for $Re_{b}=5600$ at $\phi=10^{-3}$, which is due to the collapse of fluid turbulence. In Fig. 6(c) (b and c), $\langle F_{p}\rangle_{s}$ in simulation $B_{2}$ is almost 65% lower compared to $B_{1}$. Thus, it is observed that a significant decrease in the feedback force (scaled with fluid bulk velocity ($\bar{u}$) and channel width ($2\delta$)) occurs with an increase in the Reynolds number.

The mean filtered kinetic energy equation for the mean flow is written as,

[TABLE]

In Eqn. 20, on the left-hand side, the first two terms are the energy flux terms due to molecular viscosity, eddy viscosity, and Reynolds stress. The third term is the energy utilized for the turbulent production, the fourth term is the dissipation of turbulent mean kinetic energy, the fifth term is the dissipation caused by the particle feedback, and the last term is the energy input due to pressure work. The spatial average of Eqn. 20 in the wall-normal direction makes the flux term zero due to no-slip condition at the walls, and the equation becomes,

[TABLE]

The terms in the Eqn. 21 are plotted in Fig. 7(c) for different channels and volume fractions for $Re_{b}=5600$. All the terms in Fig. 7(c) are normalized by respective $\bar{u}^{3}/2\delta$. The turbulent production term is represented by $P$, mean viscous dissipation by $\epsilon_{m}$, pressure work by $\Pi$, and dissipation due to particles by $D$. The wall-normal distance ($y$) is scaled by channel width ($h$). For the unladen case in Fig. 7(c) (a), all predicted terms for all the channels with different dimensions overlap as expected. However, for a particle volume fraction of $2\times 10^{-4}$ in Fig. 7(c) (b), the deviation is observed in all the terms. The peak of turbulent production terms for the case $C_{1}$ is lowered by nearly 30%, viscous dissipation in the near-wall is lower by 8%, and pressure work in the channel center is lower by 10% compared to the case $A_{1}$. However, $\langle D\rangle_{s}$ is higher by nearly 400% for the case $C_{1}$ compared to the case $A_{1}$. For $\phi=5\times 10^{-4}$ in Fig. 7(c) (c), the turbulence collapse in case of $C_{1}$, and the turbulent production term becomes zero. The viscous dissipation in the near-wall is lower by approximately 40%, and the pressure work is also reduced by 30%. The dissipation due to particles is higher by nearly 400% for the case $C_{1}$ than the case $A_{1}$.

In Fig. 8(c), the terms from Eqn. 21 are plotted for cases $B_{1}$ and $B_{2}$, and for different volume fractions to analyze the effect of Reynolds number. Fig. 8(c) (a) shows that for unladen flow, the peak values of turbulent production and viscous dissipation (in the near-wall region) are higher by nearly 50% for the case $B_{2}$ ($Re_{b}=13750$) compared to the case $B_{1}$ ($Re_{b}=5600$). However, $\langle\Pi\rangle_{s}$ is approximately 25% less for $Re_{b}=13750$ than $Re_{b}=5600$. This decrease in $\langle\Pi\rangle_{s}$ is compensated by the decrease in $\langle\epsilon_{m}\rangle_{s}$, while there is almost no variation in the $\langle P\rangle_{s}$. The peak value of turbulent production is much closer to the wall for a higher Reynolds number. In Fig. 8(c) (b) for $\phi=5\times 10^{-4}$, the peak value of production term and viscous dissipation in the near-wall region is almost 50% for $Re_{b}=5600$ compared to $Re_{b}=13750$. However, $\langle\Pi\rangle_{s}$ is nearly 20% and $\langle D\rangle_{s}$ is nearly 80% lower for $Re_{b}=13750$ compared to $Re_{b}=5600$. In Fig. 8(c) (c) for $\phi=10^{-3}$, the mean production term almost becomes zero, and the mean viscous dissipation in the near-wall region is 50% lower for the case $B_{1}$ compared to the case $B_{2}$. The scaled $\langle\Pi\rangle_{s}$ is nearly 5%, while $\langle D\rangle_{s}$ is nearly one order of magnitude lower for $B_{2}$ compared to $B_{1}$. From the above analysis, it is clear that the variation in particle feedback term is more significant than the others before the critical loading for any of the cases. The effect of particle feedback increases with an increase in $2\delta/d_{p}$ for fixed Reynolds numbers, and decreases with an increase in fluid Reynolds number while keeping $2\delta/d_{p}$ fixed.

III.4 Discussion

In the above section, we observed that the particle-induced drag force depends on the channel dimension and the fluid phase Reynolds number. Such an increase in drag and drag-induced dissipation is the source of different extents in turbulence modulation when $\phi$ is less than CPVL. In this section, we discuss this dependence in detail. From Eqn. 2, the scaled feedback force on the fluid due to particles can be written as,

[TABLE]

[TABLE]

[TABLE]

Here, $\phi_{l}$ is the local volume fraction of the particle, $f$ is the inertial correction factor, $Re_{p}$ is the particle Reynolds number, $Re_{b}$ is the fluid bulk Reynolds number, $\delta$ is the half-channel width, and $d_{p}$ is the particle diameter. In the case of fixed $Re_{b}$ and $d_{p}$, when $2\delta/d_{p}$ is modified by varying $\delta$, it is observed that there is no significant variation in $Re_{p}$, $\phi_{l}$, and $f$ over a range of volume fractions as shown in Fig. 9(b). Under those conditions, we can write,

[TABLE]

Thus, a change in $2\delta/d_{p}$ from 54 to 117 results in almost one order of magnitude change in the $(2\delta/d_{p})^{3}$, which is much higher than the variation of other parameters. While in the case of constant $2\delta/d_{p}$ when $Re_{b}$ is changed, there is no significant variation in $Re_{p}$, $\phi_{l}$, and $f$ over a range of volume fractions (Fig. 10(b)). Then, using Eqn. 24 we can write,

[TABLE]

Here, an increase in $Re_{b}$ from 5600 to 13750 leads to a nearly six-fold decrease in the $1/Re_{b}^{2}$ term, which is higher than the variation of other parameters. Therefore, the above formulation indicates that $2\delta/d_{p}$ and $Re_{b}$ are the dominating factors in controlling the extent of turbulence modulation compared to the other parameters.

III.5 Two point statistics: Spatial correlation

The spatial correlation of the fluid fluctuations qualitatively describes the length scale associated with the structure of the fluid phase. The spatial correlation of fluid fluctuation in the streamwise ($x$) direction is defined as,

[TABLE]

Where $r$ is the distance between the two points, and $i$ and $j$ are the components of fluctuations. The overbar denotes the ensemble averaging. The spatial correlation coefficient is defined as,

[TABLE]

$R_{ij}$ is plotted for streamwise and spanwise directions in the near-wall ($y^{+}=15$) regime for a range of particle volume fractions and shown in Fig. 11(f). In Fig. 11(f) (a and b), the correlation coefficient is plotted for $\phi=0$, $1.0\times 10^{-3}$, and $1.3\times 10^{-3}$ for the case $A_{1}$. For the streamwise direction, the correlation coefficient of the streamwise component ($R_{xx}$) decays slower than the $R_{yy}$ and $R_{zz}$ as shown Fig. 11(f) (a). The decay of the correlation coefficient significantly decreases with an increase in the volume fraction from $\phi=0$ to $1.3\times 10^{-3}$ for each of the components. However, the variation in decay rate is more pronounced for the streamwise component. Such an observation indicates that the turbulent structures become lengthier with an increase in volume fraction. A decrease in the correlation coefficient for all the components is observed when $\phi=1.5\times 10^{-3}$ (for the case $A_{1}$) as the turbulence has collapsed for this volume loading (not shown here). For the spanwise spatial correlation coefficient in Fig. 11(f) (b), no significant variation is observed in the correlation coefficient for wall-normal and spanwise components. However, an increase is observed for the streamwise component with an increase in volume loading, which signifies that the vortical structures become wider in the spanwise direction.

In Fig. 11(f) (c and d), the effect of channel dimension on the correlation coefficient is analyzed for $Re_{b}=5600$ and $\phi=5\times 10^{-4}$. In the streamwise direction for streamwise and spanwise velocity fluctuations, the correlation decay becomes slower with an increase in $2\delta/d_{p}$, Fig. 11(f) (c). Such an observation is a signature of the lengthier structure with increased channel dimensions for a fixed Reynolds number. However, in the case of correlations in the spanwise direction, the wall-normal fluctuations decay faster, and there is almost no significant modification of the dimension of turbulent streaks. The effect of Reynolds number on the correlation coefficients for the cases $B_{1}$ and $B_{2}$ are plotted in Fig. 11(f) (e and f). Interestingly, the correlation coefficient decay faster for the higher Reynolds number of $Re_{b}=13750$ compared to the $Re_{b}=5600$ for both the directions and all the components of velocity fluctuations. This suggests that the vortical structures are smaller in dimensions at a high Reynolds number. In this context, it is worth mentioning that for $2\delta/d_{p}=20$, $\phi=2.36\times 10^{-3}$ and $Re_{p}<50$, Yu et al. (2021a) have reported less turbulence attenuation for higher Reynolds number of $Re_{b}=12000$ than the $Re_{b}=5746$. They have commented that at a high Reynolds number, the large-scale vortices, although smaller in size, are stronger than those of a low Reynolds number. Zhou et al. (2020) have also reported that the vortical structures become weaker and fewer with an increase in mass loading. The authors also found that the spacing between the streaks increases with the addition of particles.

In Figs. 12(b)-13(b), the contours are plotted for instantaneous streamwise velocity fluctuations ($\widetilde{u}^{\prime}_{x}$) in the near-wall region ($y^{+}=15$) to depict the size and strength of the vortical structures. The streamwise fluctuations are scaled with respective unladen friction velocity ($u_{\tau}$). In Fig. 12(b) (a), the low and high-speed streaks are shown in the $x-z$ plane for the case $A_{1}$ at lower volume fraction ($4\times 10^{-4}$). The streamwise and spanwise dimensions are provided in viscous units. It is observed that the length of the high-speed streaks (red regime) is around 2000 viscous units, whereas low-speed streaks (blue regime) are extended up to complete channel length. For the same channel in Fig. 12(b) (b), the streaks become fewer and wider as the volume fraction is increased to $10^{-3}$. The high-speed streaks become more elongated and extend up to 3500 viscous units, while there is no significant change in the lengths of low-speed streaks. The effect of $2\delta/d_{p}$ on the high-low speed streaks for a volume fraction of $\phi=4\times 10^{-4}$ is shown in Figs. 12(b) (a) and 13(b). It is observed that the high-speed streaks become wider and fewer as the $2\delta/d_{p}$ is increased while keeping $\phi$ constant. The length of high-speed streaks extends up to 3500 viscous units from case $A_{1}$ to $C_{1}$.

In Fig. 14(d), the normalized contours for instantaneous streamwise fluctuations are shown for cases $B_{1}$ and $B_{2}$. The contours are plotted in the $x-z$ plane at the near-wall region ($y^{+}=15$) for unladen flow and a volume fraction of $5\times 10^{-4}$ while channel dimensions are kept constant at $2\delta/d_{p}=81$. In Fig. 14(d) (a and b) for the case $B_{1}$, it is observed that the streaks become fewer, lengthier, and fatter with an increase in particle volume loading. However, the streaks are smaller, thinner, and closely spaced for $Re_{b}=13750$ when compared to similar cases in Fig. 14(d) (a and c) and Fig. 14(d) (b and d). It is interesting to note that there is almost no change in the streaks spacing and lengths in case of higher Reynolds number, Fig. 14(d) (c and d). While streaks become fewer, fatter and the spacing becomes more for lower $Re_{b}=5600$ when volume fraction is increased from 0 to $5\times 10^{-4}$.

The coherent structures responsible for the turbulent production are analyzed for different $2\delta/d_{p}$, Reynolds numbers and volume fractions. Different methods have been reported in the literature to analyze the turbulent structures Chong, Perry, and Cantwell (1990); Jeong and Hussain (1995); Chakraborty, Balachandar, and Adrian (2005); Dritselis and Vlachos (2011); Yu et al. (2021c). In the present work, we report the contours of $-\lambda_{2}$ which is the second largest eigenvalue of $S^{2}+\Omega^{2}$. Here, $S$ is the strain rate tensor, and $\Omega$ is the rotational tensor of the velocity gradient. This method was proposed and discussed in detail by Jeong and Hussain (1995). This method neglects the unsteady straining which can cause minimum pressure without the presence of vortical structures and viscous effects that may eliminate minimum pressure in vortical flowsJeong and Hussain (1995). The isocontours for $-\lambda_{2}^{+}=0.0025$ shows the effect of volume fraction, system size, and Reynolds numbers, Fig. 15(c) and Fig. 16(d). The (+) symbol indicates that the values are reported in viscous units. The contours are colored according to the fluid mean velocity ($\bar{u}$) normalized with respective unladen friction velocity ($u_{\tau}$). In Fig. 15(c) (a) for case $A_{1}$, the $-\lambda_{2}^{+}=0.0025$ structures are shown for lower volume fraction of $4\times 10^{-4}$. Here, the length and width of the structures are approximately 200 and 50 viscous units, respectively. However, the $-\lambda_{2}^{+}$ structures vanish as the volume fraction reaches near CPVL (not shown here). In Fig. 15(c) (a, b and c), It is observed that the $-\lambda_{2}^{+}$ structures become fewer as channel dimension is changed from case $A_{1}$ to case $C_{1}$ while volume fraction and Reynolds number are kept constant. In Fig. 16(d), $-\lambda_{2}^{+}$ contours are plotted for both the unladen and particle-laden flows for cases $B_{1}$ and $B_{2}$. The $-\lambda_{2}^{+}=0.0025$ contours become smaller and thinner as the Reynolds number is increased from 5600 to 13750 for unladen flows. When particle volume fraction is $5\times 10^{-4}$, the number of $-\lambda_{2}^{+}$ structures becomes significantly smaller for low Reynolds number.

III.6 Particle statistics

In the previous sections, the effect of $\phi$, $2\delta/d_{p}$, and $Re_{b}$ are studied on the fluid phase properties. It is observed that the turbulence attenuation increases with an increase in $\phi$ up to a volume fraction, and a sudden turbulence collapse is observed at CPVL. However, with an increase in the $2\delta/d_{p}$ ratio, a higher extent of turbulence attenuation is observed for the fixed $\phi$ and Reynolds number. But, for the same volume fraction and fixed channel dimension, the extent of attenuation is low at a higher Reynolds number. This section discusses the modification of particle phase properties for different $\phi$, $2\delta/d_{p}$, and fluid phase Reynolds number.

The effect of the channel dimension ($2\delta/d_{p}$) on the particle phase properties at $Re_{b}=5600$ is shown in Fig. 17(f). Mean velocities and second moments of the fluctuating velocities are plotted for cases $A_{1}$, $B_{1}$, and $C_{1}$ for two different volume fractions of $\phi=2\times 10^{-4}$ and $5\times 10^{-4}$. The normalized particle mean velocities are plotted along with unladen fluid mean velocity, shown in Fig. 17(f) (a). The particle mean velocity profile is almost flat across the channel width for the case $A_{1}$. However, the particle mean velocity is lower in the near-wall region for larger channel sizes of case $B_{1}$ and $C_{1}$. For case $C_{1}$, the particle mean velocity is less in the near-wall region for $\phi=2\times 10^{-4}$ than the other volume fraction cases. At the channel center, the particle mean velocity is higher for both the volume fractions for the case $C_{1}$ as the considered Stokes number is less here than the other two system sizes, Table 1. With an increase in solid volume fraction for a fixed channel dimension, near-wall particle velocity increases due to the transfer of momentum from the channel center to the wall through particle-particle interaction. With the increase in channel dimension, Stokes number of the particle decreases, leading to a decrease in the fluid-particle relative velocity near the wall.

The normalized particle fluctuations are plotted in Fig. 17(f) (b - d). The streamwise particle fluctuations shown in Fig. 17(f) (b) are higher in the near-wall region and decrease away from the wall. The streamwise particle fluctuations decrease as the volume fraction is increased from the $\phi=2\times 10^{-4}$ to $5\times 10^{-4}$ for a fixed system size. The wall-normal and spanwise particle fluctuations increase with an increase in the volume fraction due to an increase in the collision frequency for a fixed system size Yamamoto et al. (2001), Fig. 17(f) (c and d). It is observed that the streamwise particle fluctuation increase with an increase in system size. Similar behavior is observed for the wall-normal and spanwise particle fluctuations. A more significant increase in wall-normal and spanwise particle fluctuations is observed for the case $C_{1}$ as the volume fraction is increased from $\phi=2\times 10^{-4}$ to $5\times 10^{-4}$. The dimensional particle streamwise fluctuations are plotted in Fig. 17(f) (f) for similar cases as in Fig. 17(f) (b). Here, it is observed that the dimensional particle streamwise fluctuations decrease with an increase in system size for the same volume fraction. This is similar to the behavior observed for fluid bulk velocity which decreases with an increase in system size for a constant Reynolds number. However, it is interesting to note that the particle fluctuations increase relative to the fluid bulk velocity with an increase in channel dimensions, Fig. 17(f) (b).

Next, we have analyzed the mechanisms of the generation of particle streamwise fluctuations which originate from three sources Goswami and Kumaran (2010, 2011). The first is due to the force exerted by the fluid fluctuations in the streamwise direction, which is proportional to the $D_{xx}\tau_{p}$ where $D_{xx}$ is the velocity space diffusion coefficient due to fluid fluctuations. Velocity space diffusion coefficient ($D_{ij}$) for the particle in turbulent field is defined as,

[TABLE]

Here, $R_{ij}$ is the Eulerian time correlation tensor of the fluid fluctuations, and $\int\limits_{0}^{\infty}dt^{\prime}R_{ij}$ is the integral time scale of the fluid phase ($\tau_{I}$). Eqn. 29 can be written as,

[TABLE]

In the above equation, $\tau_{I}$ may be replaced with the time scale defined as $2\delta/\bar{u}(=\tau_{f})$. Therefore,

[TABLE]

Normalizing the above relation with fluid average velocity($\bar{u}$), we get,

[TABLE]

Here, $St$ is the particle Stokes number based on the fluid integral time scale ($2\delta/\bar{u}$). The second contribution to the streamwise particle velocity fluctuations is due to particle migration across the streamlines due to wall-normal fluid fluctuations, which is proportional to the $(\tau_{p}D_{yy})St_{\gamma}^{2}$. The third source is due to collision induced by the mean velocity gradient of the particle phase, which is proportional to the $\phi_{l}(\dot{\gamma}d_{p})^{2}St_{\gamma}^{3}$. Here, $\phi_{l}$ is the local particle volume fraction, $\dot{\gamma}$ is the mean velocity gradient of particle phase, and $St_{\gamma}$ is the particle Stokes number based on the mean strain rate as $St_{\gamma}=\tau_{p}\dot{\gamma}$. Simplifying, $St_{\gamma}=\tau_{p}(dv_{p}/dy)=\tau_{p}(dv_{p}^{*}/dy^{*})(\bar{u}/\delta)=2\times St(dv_{p}^{*}/dy^{*})$, where $v_{p}^{*}$ is the normalized particle velocity, and $y^{*}$ is the normalized wall-normal distance. In the near-wall region, the streamwise particle fluctuations generated due to fluid streamwise fluctuations, which are inversely proportional to $St$, are low as the $St$ is large for cases $A_{1}$, $B_{1}$ and $C_{1}$. The contribution due to particle collisions is also low due to smaller local volume fraction ($\phi_{l}$). The significant contribution is due to the wall-normal migration of the particle, which is proportional to $St_{\gamma}^{2}$. In the channel center, the contribution due to wall-normal migration decreases due to low $St_{\gamma}$. However, it remains dominant compared to the other two sources. The streamwise particle fluctuations increase with an increase in system size (shown in Fig. 17(f) b) as $St_{\gamma}$ increases from case $A_{1}$ to $C_{1}$.

The wall-normal particle fluctuations are also generated due to three mechanisms. First, due to wall-normal fluid fluctuations, which is proportional to $(\tau_{p}D_{yy})$. The second mechanism is the collision of the particles having different mean velocities, which is proportional to the $\phi_{l}(\dot{\gamma}d_{p})^{2}St_{\gamma}$. The third contribution is due to the collision induced by the streamwise particle fluctuations, which is proportional to the $\phi_{l}(T^{1/2}_{xx}\tau_{p}/d_{p})T_{xx}$. Here, $T_{xx}=\overline{v^{\prime 2}_{x}}$ is the mean square streamwise particle fluctuations. Non-dimensionalising the $\phi_{l}(T^{1/2}_{xx}\tau_{p}/d_{p})T_{xx}$ with the fluid bulk velocity and channel width ($2\delta$), it becomes $\phi_{l}(T^{*3/2}_{xx}\times St\times 2\delta/d_{p})$. In the present study, the particle Stokes number is high, thus, the contribution from $\tau_{p}D_{yy}$ is low. The second contribution, $\phi_{l}(\dot{\gamma}d_{p})^{2}St_{\gamma}$, is also not significant as local volume fraction ($\phi_{l}$) and $(\dot{\gamma}d_{p})$ are low. The significant contribution is due to the third mechanism for $A_{1}$, $B_{1}$, and $C_{1}$. This dominant term increases (due to an increase in $T_{xx}$ and $2\delta/d_{p}$) with the increase in system sizes from $A_{1}$ to $C_{1}$ which is also observed in Fig. 17(f) (c).

The normalized particle concentration is plotted in Fig. 17(f) (e). The concentration profiles are almost flat for lower system sizes of $2\delta/d_{p}=54$ and 81 for the range of volume fractions of $\phi=2\times 10^{-4}$ and $5\times 10^{-4}$ considered in this study. However, the particle concentration is larger in the near-wall and channel center region for case $C_{1}$. This might be due to the decrease in effective particle inertia with the reduction of Stokes number (Table 1).

The effect of the Reynolds number on the particle properties is examined in Fig. 18(f) for cases $B_{1}$ and $B_{2}$. The particle properties are plotted for two volume fractions of $\phi=5\times 10^{-4}$ and $10^{-3}$. The particle mean velocity predicted for both the volume fractions are almost similar as the Stokes number is also the same as shown in Fig. 18(f) (a). The particle mean velocity is almost 10% lower in the near-wall region compared to the channel center. In case of streamwise particle fluctuations, the fluctuations are higher for $Re_{b}=13750$ compared to the $Re_{b}=5600$ for the same volume fractions near the wall. However, away from the wall, the streamwise fluctuations are lower for case $B_{2}$ ($Re_{b}=13750$) than that in case of $B_{1}$ ($Re_{b}=5600$). For particle wall-normal and spanwise fluctuations, it is observed that the fluctuations are lower for the case $B_{2}$ compared to the case $B_{1}$. The normalized particle concentration is almost flat and the same for all the considered volume fractions as shown in Fig. 18(f) (e). The dimensional particle streamwise fluctuations are plotted in Fig. 18(f) (f) for similar cases as in Fig. 18(f) (b). Here, it is observed that the dimensional particle streamwise fluctuations increase with an increase in Reynolds number for the same volume fraction. The fluid bulk velocity also increases with an increase in Reynolds number for constant system size; however, the particle fluctuations scaled with the fluid bulk velocity decrease with an increase in Reynolds number, Fig. 18(f) (b).

When volume fraction is $5\times 10^{-4}$, the component $(\tau_{p}D_{yy})St_{\gamma}^{2}$ increases in the near-wall region for the case of $B_{2}$ due to high $St_{\gamma}$. Thus, the streamwise particle fluctuations caused by wall-normal migration are higher in the near-wall region for the case $B_{2}$ compared to $B_{1}$ for $\phi=5\times 10^{-4}$ (Fig. 18(f) b). $(\tau_{p}D_{yy})St_{\gamma}^{2}$ decreases in channel center due to reduction in $\dot{\gamma}$. In the channel center, $(\tau_{p}D_{yy})St_{\gamma}^{2}$ provides almost simialr values for the cases $B_{1}$ and $B_{2}$ as $St_{\gamma}(=2\times St(dv_{p}^{*}/dy^{*}))$ is also similar for both the cases. However, it is not observed in Fig. 18(f) (b), and the particle streamwise fluctuations are lower for $Re_{b}=13750$ compared to the $Re_{b}=5600$. The reason behind this observation is not clear at this moment. For wall-normal particle fluctuations, $\phi_{l}(T^{*3/2}_{xx}\times St\times 2\delta/d_{p})$ is the dominant source, and it is lower for the $B_{2}$ as the $T^{*}_{xx}$ is lower compared to $B_{1}$ (except very close to the wall) as shown in Fig. 18(f) (b). The above discussion explains the statistics of the particle phase shown in Fig. 17(f) - 18(f).

IV Conclusion

In the present work, the simulations are performed for vertical channel flow for different volume fractions, system sizes, and Reynolds numbers to critically assess those parameters’ effects on turbulence modulation. The fluid phase is described in the Eulerian framework using the large eddy simulation. A dynamic one-equation model has been adopted to model the subgrid-scale stress. The individual particles are tracked using a Lagrangian technique with point-particle assumption. It is observed that the fluid fluctuations decrease with an increase in volume fraction up to a critical particle volume loading at which turbulence collapses. The extent of turbulence attenuation increases with an increase in system size ($2\delta/d_{p}$) for the same volume fraction while keeping the Reynolds number fixed. But, for the same volume fraction and fixed channel dimension, the extent of attenuation is low at a higher Reynolds number.

Different terms from the momentum and energy balance equations are plotted, and it is concluded that the variation of the particle feedback term is the dominant term when the volume fraction is lower than the critical volume fraction. It is observed that the dissipation due to the particle significantly increases with an increase in system size for the same volume fraction and fixed Reynolds number. However, the dissipation due to particle decreases significantly when the Reynolds number is increased while keeping channel dimensions same, which reduces the extent of turbulence attenuation. From the analysis of the variation of scaled feedback force, it is found that the system size and bulk Reynolds number are the two important parameters that significantly affect the turbulence modulation compared to the others.

The spatial correlation of fluid fluctuations is plotted in the streamwise and spanwise directions, and it is observed that the correlation coefficient for the streamwise component decays slower with an increase in system size, which indicates that the streamwise turbulent structures become elongated along the channel length. The isocontours of the streamwise fluctuations also depict a similar observation. The isocontours show that the high-speed streaks become elongated and fewer with an increase in system size for a fixed Reynolds number. This may be because the presence of particles affects the breakdown of coherent structures with an increase in system size. The correlation coefficient decays faster for a higher Reynolds number than a lower one. The isocontours of the streamwise fluctuations show that the streaks are smaller, thinner, and closely packed for a higher Reynolds number.

The particle fluctuations, scaled with fluid bulk velocities, are analyzed for the effect of system size and Reynolds number. For a constant system size and Reynolds number, the streamwise particle fluctuations decrease, while wall-normal and spanwise particle fluctuations increase with an increase in particle volume fraction due to the increase in particle collision frequency. It is observed that the particle fluctuations increase with an increase in system size for a constant Reynolds number and volume fraction. However, the particle fluctuations, except for the streamwise particle fluctuations in the near-wall, decrease with an increase in Reynolds number for fixed channel dimension and volume fraction. Qualitative analysis for the source of particle streamwise and wall-normal fluctuations is discussed. It is observed that the streamwise particle fluctuations generated due to wall-normal migration are dominant compared to the other two sources (force due to streamwise fluid fluctuations and collision-induced fluctuations). The source due to particle migration in the wall-normal direction increases with an increase in system size. For the wall-normal fluctuations, the source due to collision induced by streamwise particle fluctuations is dominant in the near-wall and channel center. The present analysis highlights the multiscale nature of particle-laden flows in the sense that an increase in Stokes number of particles universally demonstrates the turbulence modulation, but the system size is to be considered. Therefore, this study will be a guideline for designing industrial instruments while using them for gas-solid operations.

V Appendix

The unladen fluid statistics for $Re_{b}=5600$ is presented in Fig. 19(d). The simulations are verfied against inHouse DNS Muramulla et al. (2020) and DNS data of Moser et al. (1999). The fluid phase velocities are normalized with unladen friction velocity ($u_{\tau}$), and the wall-normal distance is normalized with $u_{\tau}$ and half-channel width ($\delta$). The streamwise mean velocity and fluctuations plotted in the wall-normal direction, are in good agreement with both the DNS datas as shown in Fig. 19(d) (a and b). The cross-stream stress and wall-normal fluctuations are underpredicted by LES within 15% and 20% deviations compared to DNS datas, respectively.

In Fig. 20(d), the unladen fluid statistics for $Re_{b}=13750$ is validated against the DNS data of Moser et al. (1999) at $Re_{\tau}=395$, experimental data of Li et al. (2012) at $Re_{\tau}=430$, and LES data of Salmanzadeh et al. (2010) at $Re_{\tau}=350$. The mean velocity and fluid fluctuations predicted by LES model are in good agreement. The peak value of fluid fluctuations are within 15% deviation compared to DNS data of Moser et al. (1999).

Acknowledgement

We would like to thank the Science and Engineering Research Board (SERB), Department of Science and Technology (DST), Government of India, for their financial support.