Analysis of Mobility Effects in Particle-Gas Flows by Particle-Resolved LBM-DEM Simulations

Particle-resolved direct numerical simulations are performed to simulate the flow through particle assemblies that are either static or freely moving to demonstrate the influence of particle mobility. To obtain a comprehensive understanding for this influence essential parameters such as the Reynolds number, solids volume fraction, particle-fluid density ratio, collision parameters and particle shape are varied. The influence of particle mobility is assessed by evaluating the particle-fluid forces, the particle ensemble structure and particle velocities. It is found that the ability of existing correlations for static particle systems to predict drag and lift forces correctly in dynamic particle-gas flows is limited and that drift forces perpendicular to the drag force play an important role.


Introduction
Fluidized systems occur in a wide range of industrial application such as drying, coating, granulation, gasification, or combustion. In gas-solid fluidization the fluidized bed is characterized as a gas-solid mixture that behaves like a fluid due to the introduction of an upward gas flow that counteracts the solid particles' gravity. Numerical simulations of fluidization processes play an increasingly important role in the design and optimization thereof and can be divided into three different classes: microscopic particle-resolved direct numerical simulations (PR-DNS), mesoscopic Euler-Lagrange DEM-CFD simulations and macroscopic twofluid model (TFM) simulations [1]. PR-DNS resolve the fluid flow around individual Lagrangian particles and allow a direct calculation of the particle-fluid forces acting upon the particles whose motion can be predicted by the discrete element method (DEM). Unresolved DEM-CFD simulations also consider the motion of individual particles, but rely on using a computational fluid dynamics (CFD) method to solve the volume-averaged Navier-Stokes equations in which the assumed fluid cell size is usually larger than the particle diameter. In TFM simulations both, the fluid and the solid phase, are solved using volume-averaged equations and the motion of single particles is not tracked.
While DEM-CFD and TFM simulations are generally computationally cheaper than PR-DNS, their drawback is that they rely on constitutive models for the interphase drag forces, which have a large influence on the simulation's accuracy. In the past many different drag force correlations have been proposed for the flow through random monodisperse assemblies of spheres for varying solids volume fractions and Reynolds numbers based on PR-DNS where the particles are assumed to be rigid [2][3][4][5][6][7]. Only recently it has been questioned whether these static correlations are also accurate for freely moving particles, which are ubiquitous in industrial applications. By performing PR-DNS Tang et al. [8] and Huang et al. [9] have shown that particle velocity fluctuations increase the interphase drag at intermediate Reynolds numbers and proposed drag correlations based on the systems' granular temperature, which is a measure for the ensemble-averaged particle velocity fluctuations. Rubinstein et al. [10], on the other hand, concluded from particleresolved lattice Boltzmann (LB) simulations that in low Reynolds number flows drag forces are lowered by the ability of mobile particles to adjust to the surrounding fluid flow and proposed a Stokes number dependent drag correlation for the Stokes regime that takes into account the solid-fluid density ratio. Tavanashad et al. [11] also used LB simulations to study the effect of particle velocity fluctuations for varying density ratios at Re = 20 and found that the density ratio of mobile particles has a negligible effect on the mean drag force. From these studies it can be concluded that the mobility of spherical particles has a different impact on the drag forces depending on the flow regime. However, a generalized theory or model that explains the observed phenomena over a wide range of flow regimes is still missing.
While many previous studies focused on interphase drag forces, only few have addressed lift forces in fluidization simulations. Esteghamatian et al. [12] compared PR-DNS and DEM-CFD simulations of bubbling gas-solid fluidized beds and found that available analytical models for lift force contributions such as the Saffman force and the Magnus force can neither predict the high magnitude nor the direction of the lift force computed by DNS correctly. Moreover, Esteghamatian et al. [12] found that the particles' transverse slip velocities and fluctuations in the transverse force component are both underpredicted by DEM-CFD simulations. A contrary observation was made by Liu and Wachem [13], who compared fluidization experiments in the bubbling regime with corresponding DEM-CFD simulations and reported an overestimation of time-averaged horizontal particle velocities and horizontal velocity fluctuations in their simulations. Thus, the correct implementation of lateral forces in DEM-CFD simulations and their impact on the simulation accuracy are still an open question.
Until now, lift force models based on DNS results have mainly been developed in the context of non-spherical particles. However, while drag and lift force correlations with variable angle of attack have been proposed for isolated non-spherical particles [14][15][16][17], there is generally a lack of correlations for dense non-spherical particle systems. A rigorous drag and lift force correlation for spherocylinders with aspect ratio 4 that is valid for variable Reynolds numbers and solids volume fractions has been recently developed from particle-resolved LB simulations of flow through static spherocylinder assemblies by Sanjeevi and Padding [18]. The authors of the latter study highlighted the importance of lift force, since their magnitude can especially at higher Reynolds numbers attain values as high as the drag force [18]. The influence of mobility effects, i.e., effects related to the ability of particles to change their position and velocity in a flow as a result of particle-particle and particle-fluid interactions, on the accuracy of static interphase drag or lift correlations for non-spherical particles has not been thoroughly addressed yet. Therefore, it remains unknown if the recent findings for mobile spherical particles are also transferable to particles of complex shapes.
In the present study we aim at extending the understanding for mobility effects in dense gas-solid fluidized systems by performing particle-resolved lattice Boltzmann (PR-LB) simulations of flow through both, static and freely moving assemblies of particles, which are coupled with discrete element modelling to account for particle collisions that frequently occur when particles are allowed to freely move. PR-LB simulations as a subgroup of PR-DNS are in this case a very promising approach, because they allow for analyzing the directly computed particle-fluid forces of individual particles that are not accessible in unresolved DEM-CFD simulations. Note that in the following to the performed PR-LB simulations will be referred as PR-DNS. The following parameters that are expected to have an influence on the mobility effects are varied in the present study: Reynolds number, solids volume fraction, particlefluid density ratio, particle shape and collision behavior, i.e., restitution coefficient and friction coefficient. Especially the effects of particle shape and collision parameters have not been addressed in DNS studies on mobility effects yet, but could reveal and explain further limitations of static correlations. This work is divided as follows. In Sect. 2 the applied coupling between the lattice Boltzmann method and the discrete element method are briefly outlined and the performed simulations are explained. In Sect. 3 the observed mobility effects are addressed by considering three aspects: deviation from static drag as well as lift force correlations, particle ensemble structure (relative position of particles), and particle velocity correlations. Finally, in Sect. 4 the conclusion, that motivates further investigations of mobility effects in future works, is drawn.

LBM-DEM Coupling
The particle-resolved direct numerical simulations of the fluidized systems investigated in this work were carried out using a coupling of a lattice Boltzmann flow solver and the discrete element method. The coupling of used lattice Boltzmann flow solver and DEM solver has been previously validated for suspension flows with moving particles in Rosemann et al. [19] and Kravets et al. [20]. The LB algorithm is based on the collision of velocity distribution functions in fluid nodes separated by the grid spacing Δx and the consecutive streaming of post-collision distribution function to neighboring fluid nodes. This fundamental process can be written as where the distribution functions f i are streamed in one time step Δt from the node at position r along the discrete lattice velocities e i to the neighboring nodes at position r + e i Δt by taking into account the collision operator Ω(f) and external forces F i . Here, a three-dimensional D3Q19 lattice model is used that features N = 19 discrete lattice velocities and the multiple relaxation time (MRT) collision operator presented in d'Humièrs et al. [21], which allows using different relaxation times for different moments of the Navier-Stokes equation besides the viscosity-linked relaxation time t v to enhance the stability of the method. The no-slip boundary condition for particles is established by the linearly interpolated bounce-back method by Bouzidi et al. [22] that computes the ''bounce-back'' of distribution functions at the particle surface during a streaming step. The particle-fluid force F p-f and resulting torque T p-f are computed using the momentum exchange method by Wen et al. [23]. Solid nodes that are turned into fluid nodes in the wake of moving particles are initialized with the simple equilibrium refilling method which performs equally well as more sophisticated refilling methods in dynamic particulate flows [19]. Further details about the LBM implementation can be found in Rosemann et al. [19]. The particle trajectories are calculated by the discrete element method that tracks the particle positions x p , translational velocities u p and rotational velocities _ w p according to Newton's equations of motion given by and where the particle-particle collision force F c , the particle-fluid force F p-f and the corresponding torques T p-f are acting upon the particles with mass m p and inertia tensor I p . The collision force includes a normal and a tangential component. The normal component is modelled as a linear spring damper model which depends on the normal spring stiffness and the normal damping coefficient. These two coefficients are calculated based on the particle mass, collision duration and restitution coefficient e as outlined in Schäfer et al. [24]. A collision duration is chosen that is 25 times the DEM integration time step Δt DEM , which is necessary to resolve particle collision with sufficient accuracy [25]. For the tangential contact force, a spring model with a tangential spring stiffness that is 0.9 times the normal spring stiffness is used which is additionally constrained by the Coulomb friction based on a friction coefficient m. Rolling friction is neglected in the collision model.

Simulation Setup
To investigate mobility effects in gas solid-flows the flow through particle systems is considered as depicted in Fig. 1, where particles are either assumed to be static with zero velocity or able to freely move due to the forces and moments acting on them. All of the analyzed systems include N p = 100 particles in a cubic domain with length L and fully periodic boundaries. The solids volume fraction j is adjusted by the domain size to be either 0.2 or 0.4. The number of particles in the present work is comparable with the number of particles used in the above mentioned DNS studies investigating particle-fluid forces in assemblies of spherical [2][3][4][5][6][7][8][9][10][11] and non-spherical [18] particles, where it varies between 27 and 200 depending on the study. In this study two different particle shapes were considered: spheres and spherocylinders with aspect ratio 4 that have the same volume as the considered spheres. These spherocylinders were also used in the previous study of Sanjeevi and Padding [18], but instead of an exact spherocylinder shape representation the multi-sphere approach is used with which the spherocylinders are modelled as 11 aligned and equispaced spheres as can be seen in Fig. 1. The simulations start with a random configuration being created by two consecutive steps: initialization at overlap-free position and randomization of particle positions. These two steps are carried out for the two investigated particle shapes as follows.
The spheres are first inserted one after another at random positions in the cubic domain that do not lead to overlaps with already inserted particles. When all spheres are inserted a Monte-Carlo simulation follows in which the particle positions are varied by a small displacement in each time step that is only accepted if the new position does not lead to an overlap with an already inserted particle. The spherocylinders are first inserted at random overlap-free positions with random orientations in a larger domain with size L L H surrounded by walls, where the height H is larger than L such that each particle can be inserted with a minimum midpoint-to-midpoint distance of one cylinder length. A downward gravity force is then assigned to the particles to make them settle to the bottom of the enclosure. During the settling process friction-free particle-particle and particle-wall collisions are assumed using the discrete element method. Once the particles are settled, they fit into the cubic domain with size L L L. With the final particle positions of the settling process another DEM simulation with fully periodic boundaries in the desired cubic domain is initialized, where random translational and rotational velocities are assigned to the particles and fully elastic and friction-free collisions are assumed during the simulation. After a sufficiently long simulation period this results in the desired random particle configuration in the cubic domain.
To achieve different flow conditions, the superficial Reynolds number are adjusted.
Eq. (4) depends on the particle diameter D p (volumeequivalent diameter of a sphere in the case of spherocylinders), the fluid's kinematic viscosity v and on the slip-velocity between fluid and particles given by Where <u p > denotes the ensemble-averaged particle velocity and u f the fluid velocity averaged over all fluid nodes in the domain. The desired particle Reynolds numbers are established by applying an appropriate external fluid force Fi in Eq. (1). The total external fluid force is obtained by summing up this force contribution over all fluid nodes. In the case of freely moving particles, the total fluid force divided by the number of particles is imposed into opposite direction on each particle to keep the system's net momentum zero and achieve a constant time-averaged Reynolds number.
The LB simulation parameters are listed in Tab. 1. The resolution D p /Δx and relaxation times t v are chosen to be very close to the recent study of Sanjeevi and Padding [18], who performed a resolution-dependency study for the obtained drag forces acting upon static ensembles of spherocylinders under various flow conditions that are similar to the flow conditions considered in this study, i.e. similar solids volume fractions and Reynolds numbers. For reasons of consistency the same simulations parameters for spheres and spherocylinders are chosen. Since 100 particles are used in each simulation the resulting domain lengths L are 6.4D p for the solids volume fraction j = 0.2 and 5.1D p for j = 0.4. The LBM time step has been chosen such that the Mach number of the flow is kept below 0.3 during the simulation, which is a stability criterion of the underlying incompressi-ble lattice Boltzmann scheme. The DEM time step has to be sufficiently small to keep the overlap between particles during a collision below 1 % and obtain a correct collision behavior. Since the LBM time steps at lower Reynolds numbers (Re ( 1 and Re = 20) would be too large to ensure the required small particle overlaps, DEM time steps smaller than the LBM time steps have to be chosen and the ratio Δt LBM /Δt DEM is therefore larger than one for lower Reynolds numbers. This means that the particle-fluid force F p-f in Eq. (2) is assumed to be constant for Δt LBM /Δt DEM DEM time steps.
In this work three different simulation cases were considered to probe the influence of particle density and contact parameters on the correlations' ability to accurately predict the drag force for freely moving particles, which are given in Tab. 2. Case I and II assume a fully elastic and frictionfree collision behavior but differ in the particle-fluid density ratio (case I: r p /r f = 1000, case II: r p /r f = 100). The density ratio in case III is the same as in case II but in case III an inelastic collision behavior with restitution coefficient e = 0.5 and friction coefficient m = 0.3 is simulated. These values are exemplarily chosen here to demonstrate the possible effect of elasticity and friction. However, we note that the chosen values for the friction coefficient and the restitution coefficient are close to the corresponding values of biomaterials such as corn grains [26].  Table 2. Simulation parameters of the three simulation cases involving freely moving particles: particle-fluid density ratio r p / r f , friction coefficient m and restitution coefficient e.

Particle-Fluid Forces
The main contribution of particle-fluid forces in gas-solid flows is usually assumed to be the drag force which is aligned with the slip-velocity. In the present work the particle-fluid force is decomposed into the drag force and the perpendicular drift force as follows (see Fig. 1): Note that the particle-fluid force F p-f does only include the viscous contribution and the contribution from the fluctuating pressure field. The contribution of the external body force that is driving the flow is not included. Note also that especially in studies on interphase drag of gas-solid flows considering spherical particles F p-f and F drag are often regarded as identical since the drift force averaged over all particles is generally assumed to be zero and only the ensemble-averaged forces were of interest in these studies. As can be seen in Fig. 1, the drift force contains in the case of spherocylinders also the contribution from the lift force based on the particle orientation that cannot be defined for the spherical particle. The drag force vector can be calculated as Where F drag denotes the scalar projection of the particlefluid force onto the unit vector ê drag being aligned with the slip-velocity (Eq. (5)). The drift force is computed as with the unit vector Analogously, the lift force for the spherocylinders is defined as and the definition of ê lift given in Mema et al. [27] is used, which calculates the expected lift force direction based on the particle's orientation with respect to the slip-velocity. Note that the values F drag and F lift can be negative in case the drag force acts into the opposite direction of the slip-velocity or the lift force acts into the opposite direction of ê lift , respectively. In the following the particle-fluid forces obtained from the DNS F DNS are compared to the predicted force F cor given by reference correlations from literature. The predicted force is not based on the resolved flow around the particles but only on the superficial Reynolds number as defined in Eq. (4). The relative deviation between the two forces is computed for each particle as where the ensemble-averaged mean of the force magnitude is used as the denominator to avoid near zero divisions which would distort the deviation computations. The relative deviation of the force magnitude is calculated as follows: Tab. 3 lists the drag force deviations for the flow through a static ensemble of spheres and the three different moving cases (see Tab. 2). Note that the deviations are time-averaged for the simulations with freely moving particles and also for the static simulations at the highest Reynolds number (Re = 1000).
In the Stokes regime (Re ( 1) the static correlation by van der Hoef et al. [28] for rigid particle systems can accurately predict the drag force for the static simulation case but overpredicts the drag force for the moving cases, especially for low particle-gas density ratios. The Wen and Yu [29] correlation was obtained from sedimentation experiments and thereby inherently considers the particle mobility. A lower deviation is obtained with this correlation for the low density-ratio case II, but not for case I. The reason for this observation is that the density ratio has a notable influence on the interphase drag but is not included in the Wen and Yu [29] correlation. The recently developed correlation of Rubinstein et al. [10], on the other hand, explicitly takes into account the density ratio in form of the Stokes number St = Re/18 r p /r f and can reduce the maximum deviation observed in these cases. There is not a significant influence of friction and inelasticity during the collisions (case III) in the Stokes regime for the low solids volume fraction but a notable reduction of interphase drag for j = 0.4, which can be attributed to the fact that collisions happen more often in dense particle systems. The inelasticity and friction promote the clustering of particles and the drag reduction since the relative particle velocities are reduced during a collision event. The deviation of the correlation by Rubinstein et al. [10] is also significantly higher in case III since it was obtained solely from simulations with friction-free and perfectly elastic collisions.
At elevated Reynolds numbers (Re = 20 and Re = 1000) the static correlation of Tang et al. [7] is applicable for the simulations with static spheres and its modified version [8] for freely moving spheres is applicable to moving case I, but is outside the applicability range of case II and III because particle-fluid density ratios below 500 were not considered in the study of Tang et al. [8]. The modified version does not only take into account the slip Reynolds number (Eq. (4)) and solids volume fraction, but also the Reynolds number based on the granular temperature defined as with the granular temperature and its components Θ k that are derived from the velocity components of the particle velocities u p in each spatial direction k by the following sum over all particles n as follows: The granular temperature can be understood as a measure for the ensemble-averaged particle velocity fluctuations that occur when particles are allowed to continuously change their velocities due to the forces acting upon them. From their simulations Tang et al. [8] derived a correlation for Re T that is dependent on the particle-fluid density ratio and the Reynolds number. There is already a fair deviation between the DNS results and the correlation for the static simulation case which cannot only be attributed to different underlying numerical methods (lattice Boltzmann method vs. immersed boundary method) but also to the fact that Tang et al. [7] performed the simulations with a relatively small resolution (D p /Δx~12) and a grid-size effect correction. The number of particles cannot be the origin of the deviation since the number particles used in the study of Tang et al. [7] and the present study is approximately the same (108 vs. 100) and periodic boundary conditions are also used in both of the studies. In our previous work [6], we have extensively compared results of the used LB solver for flow through static assemblies of spheres for a broader range of intermediate Reynolds numbers and solids volume fraction with available references and found that our results agree on average better with DNS studies such as Bogner et al. [4] and Tenneti et al. [5] where higher resolutions than in the study of Tang et al. [7]. are inherently used. It has also to be noted that there is generally a lack of highly resolved simulation data for gas-solid flows at elevated Reynolds numbers, even for spherical particles.
In contrast to the Stokes regime, Tang et al. [8] found in PR-DNS that the interphase drag in dynamic gas-solid flows increases significantly with increasing Reynolds numbers and decreasing particle-fluid density ratios compared to the flow through a static particle system. However, the DNS results of the present study do not support that the effect is as large as proposed by the correction by Tang et al. [8]. The static correlation [7] predicts the drag force equally well and in most systems even much better than the dynamic correlation [8] for the moving case I. This is even  Tang et al. [8] j, Re, Re T (r p /r f , Re) more visible for the moving case II which is, however, outside the correlation's applicability range due to the low density ratio that has a very high impact on the correlation. Since the deviations between our DNS results and correlation [7] are almost identical for moving case I and II, the respective ensemble-and time-averaged drag forces are not significantly different despite the different density ratios. These finding regarding the rather low influence of the density ratio on the ensemble-averaged drag forces acting on spheres in dynamic systems are consistent with the recent findings of Tavanashad et al. [11,30], who performed similar PR-DNS for Reynolds numbers between 20 and 100 and found only significant deviations for liquid-solid systems (r p /r f 10). Regarding the moving case III, the DNS results of the present study further suggest that elasticity and friction significantly increase the drag force at intermediate and higher Reynolds numbers, which is a contrary observation to the Stokes regime where interphase drag is slightly reduced by the introduction of this collision behavior. A possible explanation for this effect could be that particles experience a sudden loss of kinetic energy during an inelastic collision with friction, which increases the slip velocity with respect to the surrounding fluid and thereby the drag force.
Tab. 4 shows the deviations between the drag and lift forces computed in the PR-DNS of spherocylinders and the respective forces obtained by the correlations of Sanjeevi and Padding [18]. These correlations do not only depend on the superficial Reynolds number (Eq. (4)) and the solids volume fraction but also the angle of attack q, which is the angle between the slip velocity u s and the longitudinal cylinder axis (see Fig. 1). First, the results for the static simulation case are discussed. Since the lift forces are very sensitive with respect to the particle orientations, the results for this case are additionally averaged over three independent particle configurations. Very low deviations with respect to the correlations are obtained regarding the values for F drag cor and F lift cor , which can be explained by the fact that the same underlying numerical method and resolution is used as in the study of Sanjeevi and Padding [18]. One difference is the number of considered particles, which is 200 in Sanjeevi and Padding [18] and 100 and in the present study. However, the good agreement suggests that the size of the system has no significant influence when studying static particle configurations. The second difference is that Sanjeevi and Padding [18] derived the correlations based on particle configurations where the particle orientation is only varied around one axis such that particles are aligned with respect to one plane parallel to the mean flow. The fact that particles have random orientations around all principal axes in the present study proves that this difference does not translate into different particle-fluid forces, which was also exemplarily shown in Sanjeevi and Padding [18]. Finally, we note that the particle shape representation by the multisphere approach does apparently not lead to different drag or lift forces compared to an exact shape representation used in Sanjeevi and Padding [18], which indicates that the number of used spheres for one spherocylinder (11 spheres) is sufficient.
When comparing the deviations for the signed forces with the force magnitudes in the static case, there is no difference between the values for the drag force (d rel (F drag cor ) = d rel (⏐F drag cor ⏐)), but a notable difference between the values for the lift force, since the correlation is significantly underpredicting the lift force magnitude. The reason for this observation is that we consider the ensemble-averaged force deviations in Tab. 4 and can be explained as follows: Especially for angles of attack close to 0°or 90°the lift force is expected to be close to zero. However, due to the deflection of the flow at neighboring particles, this is rarely the case and significant force contributions in the direction of the lift force are also obtained for these angles. These force contributions are zero on average since they have varying directions. However, when considering the force magnitudes, Table 4. Time-and ensemble averaged relative deviation of the static drag and lift force correlation of Sanjeevi and Padding [18] as defined in Eq. (11) for the simulations with spherocylinder particles being either static or freely moving (moving case II). Correlation parameters are Reynolds number Re, solids volume fraction j, and angle of attack q. Results for the static case are averaged over three independent random particle configurations. Positive values indicate an overprediction and negative values an underprediction of the DNS results by the correlation. these fluctuating force contributions do not cancel out in the ensemble average. This effect becomes less pronounced for increasing Reynolds numbers because the lift force contribution to the total particle-fluid force becomes larger with increasing Reynolds number [18]. Thus, the fluctuating force contribution in the lift force decreases.
The results for moving case II in Tab. 4 reveal that both, the lift forces and drag force, are significantly lower in the Stokes regime when the cylinder particles are allowed to move in the PR-DNS compared to the static case, which leads to an overprediction of these forces by the correlation. For higher Reynolds numbers the mean drag and lift forces do not differ significantly from the correlations. These findings are comparable with the findings for spheres in Tab. 3, but the increase in interphase drag for moving particles in the Stokes regime is more pronounced for the spherocylinders. This can be explained by the fact that particle cluster formations causing the drag reduction are more persistent when considering the elongated particles due to a higher degree of entanglement which is not possible for spherical particles. When it comes to the averaged magnitude of the forces, the deviations in the drag forces are very close (d rel (F drag cor )~d rel (⏐F drag cor ⏐)) as in the static case. The values for d rel (⏐F drag cor ⏐), however, show a different behavior compared to the static case. In the moving case the magnitude of lift forces increases with increasing Reynolds number and with increasing solids volume fraction. In the Stokes regime one can also observe that ⏐d rel (⏐F drag cor ⏐)⏐ < ⏐d rel (F drag cor )⏐ for moving particles, which is caused by a very low correlation between the predicted direction of the lift force and the real direction in the PR-DNS (F lift DNS is negative in Eq. (11) when the DNS lift force acts into the opposite direction of the predicted lift force direction). Fig. 2 shows the ratio between the ensemble-averaged magnitude of the drift force (see Fig. 1) and the ensembleaveraged magnitude of the drag force. From Fig. 2a it can be seen that the significance of the drift force for spherical particles increases with increasing solids volume fraction and is especially high in the Stokes regime. Furthermore, the ratio www.cit-journal. increases with increasing Reynolds number outside the Stokes regime and increases with decreasing particle-fluid density ratio. Note that the static case can be regarded as the limit r p /r f ¥. The ratio is only significantly influenced by inelasticity and friction during the collisions at higher Reynolds numbers (Re = 1000) and reduced by this collision behavior. Fig. 2b shows the ratio between the ensemble-averaged magnitude of the drift force and the ensemble-averaged magnitude of the drag force for simulations with static spherocylinders and freely moving spherocylinders (moving case II). This figure includes not only the PR-DNS values but also the ratio that is obtained by computing the particle drag forces and drift forces based on the correlations of Sanjeevi and Padding [18]. Note that the correlation prediction of the drift force only includes the lift force contribution and not the third contribution being perpendicular to both, lift force and drag direction, since there is no available correlation for this force contribution which purely arises from the deflection of flow in the dense particle system. In both, static and moving case, the ratio directly computed in the PR-DNS is significantly larger than the ratio predicted by the correlation and the discrepancy is more pronounced for the moving case than for the static case. The large ratios suggest that the lift force is not only underpredicted by the correlation, but also that the third force component has, especially in the moving case, a magnitude that is comparable with the lift force magnitude. When comparing the results shown in Fig. 2a and Fig. 2b, it can be seen that the ratios obtained from the PR-DNS for static case and moving case II do not differ notably. This indicates that the magnitude of the drift force relative to the drag force magnitude is not strongly shape-dependent regarding the two considered particle shapes.
The significance of the inclusion of a lift force model in unresolved DEM-CFD fluidization simulations with the spherocylinders considered in the present work has been recently studied by Mema et al. [31]. The authors found that the inclusion of lift forces has a notable effect on the flow velocity in mean flow direction (opposite direction of gravity). The results of the present study further underline the significance of drift forces in general, but suggest that lift force models based on the particle orientation alone may be insufficient to predict these force contributions correctly, especially when it comes to dynamics gas-solid flows. This assessment regarding the significance also holds for spherical particles which do not feature a specific orientation. A special attention needs to be paid to the transition from the Stokes regime to higher Reynolds numbers, because the influence of particle mobility and the magnitude of drift forces are remarkably different depending on the flow regime.
The observation that the lift force model for non-spherical particles alone is insufficient to fully describe the forces acting perpendicular to the drag force is further supported by Fig. 2c, where the ratio of the force being perpendicular to both, lift force and drag force, and the lift force is shown. It is even slightly higher than the lift force for the Stokes regime and at Re = 20, and increases significantly at Re = 1000 when particles are allowed to be mobile. This further explains the discrepancies in Fig. 2b and highlights the importance of forces that arise from the deflection of flow at neighboring particles.

Particle Ensemble Structure
To further analyze the effects of particle mobility, the underlying structure in the particle system is of interest. While static correlations usually assume a random distribution of particles, this randomness does not necessarily need to be present in a dynamic fluidized system. One tool that has been recently used to find structures in the particle system are radial distribution functions [8,11]. The radial distribution function g(⏐x p,i -x p,j ⏐) is a measure for the probability of a particle i with position x p,i to be separated by another particle j with position x p,j by the distance ⏐x p,i -x p,j ⏐. If g(⏐x p,i -x p,j ⏐) = 1, the particle separation ⏐x p,i -x p,j ⏐ in the particle system is as likely as finding the particle separation in a fully random particle distribution, since the distribution function is normalized by the domain-averaged particle density. Significantly larger values than unity for g(⏐x p,i -x p,j ⏐) indicate that particles are preferably found with separation distance ⏐x p,i -x p,j ⏐ in the particle system. When particle clustering occurs and many particles are in contact, one would expect peaks near ⏐x p,i -x p,j ⏐/D p = 1, since particles are in contact at this separation distance. Fig. 3 depicts the radial distribution functions for the moving cases I-III simulated with spheres (a-c) and the moving case II simulated with spherocylinders (d). The lower particle-fluid density ratio promotes particle clustering only in the Stokes regime and not higher Reynolds numbers outside the Stokes regime (compare Fig. 3a and 3b). The inclusion of particle elasticity and friction, on the other hand, only promotes clustering for the higher Reynolds numbers, while it is not a significant effect in the Stokes regime (compare Fig. 3b and 3c). These observations regarding the influence of the density ratio and the collision behavior depending on the flow regime are consistent with the observations made with respect to the drag forces in Sect. 3.1. Thus, the tendency for particle clustering is one possible explanation for the observed decrease (Stokes regime) or increase (higher Reynolds numbers) of the drag forces. One final observation regarding the simulations of moving spheres can be made in Fig. 3a where the simulations performed at Re = 20 have the lowest value for g(⏐x p,i -x p,j ⏐) = 1 compared to the simulations with the other two Reynolds numbers and thereby the lowest tendency for particle clustering. This can be explained by the fact that two phenomena can cause particle clustering: First, long-term particle cluster formations and void regions occur in the Stokes regime due to the particle's ability to adjust to the flow [10]. Second, at elevated Reynolds numbers strong particle velocity fluctuations introduce a more chaotic behavior of the particle system that favors frequent particle collisions. Since the particle velocity fluctuations increase with Reynolds number [8], particles can be found more often in close vicinity to each other at Re = 1000 than at Re = 20.
Interestingly, the spherocylinders do not have a peak close to ⏐x p,i -x p,j ⏐/D p = 1 but a drop in the radial distribution function (Fig. 3d). A possible explanation for this phenomenon is that the distance between the cylinder centers and not the closest distance between two cylinders is considered. It is more likely that the cylinders' touch each with outer particle surface regions than with surface regions close to the center. Especially for non-spherical particles, better tools for an analysis of the particle system structure would be desirable.
Another aspect concerning the structure in the particle system is the preferred orientation when it comes to nonspherical particles. Fig. 4 shows the probability of finding a spherocylinder with angle of attack q in the simulation with freely moving particles (moving case II). For the derivation of the static correlation of Sanjeevi and Padding [18] particle configurations with an even distribution of attack angles were used. In a dynamic system, however, the particle orientations can generally not be assumed to be evenly distributed. As the Reynolds number increases, particles are increasingly oriented perpendicular to the mean flow velocity. This observation is consistent with the results of the unresolved DEM-CFD fluidization simulations of Mema et al. [31], who found that spherocylinders are preferably aligned horizontally in the upper region of the vertical fluidization column, where particles are more mobile and no longer densely packed. This flow situation is similar to the flow situation in the simulations of the present study. Since these DEM-CFD simulations also considered gas-solid fluidization, the Reynolds numbers can be expected outside the Stokes regime and in the range of Re = 1000, where also an orientation perpendicular to the mean flow direction is obtained. The frequent occurrence of particles with attack angles between 60°and 90°supports the idea that lift forces are not strongly correlated with the angle of attack (relative to the mean flow velocity) in dynamic gas-solid flows,  . Ensemble-and time-averaged radial distribution functions for the fluidized system simulations with moving particles under different flow conditions: a) moving case I -high particle-fluid density ratio with spheres, b) moving case II -low particle-fluid density ratio with spheres, c) moving case III -low particle-fluid density with additional particle elasticity and friction with spheres, d) moving case II -low particle-fluid density ratio with spherocylinders.
because a low lift force would be predicted by a static lift force correlation for such angles [18]. However, as can be seen in Fig. 2b, the drift forces (and thereby also lift forces) have a significant influence in the simulations with freely moving particles at Re = 1000 despite high angles of attack.

Particle Velocities
Besides the particle ensemble structure, an analysis of particle velocities can also shed light on the mobility effects in dynamic gas-solid flows. Fig. 5 shows the normalized correlation between the particle velocities of two particles i and j against their separation distance ⏐x p,i -x p,j ⏐ for moving cases I-III simulated with spheres (a-c) and moving case II simulated with spherocylinders (d). Values near zero indicate that the particle velocities are uncorrelated and values significantly larger than zero indicate that the particle velocities are positively correlated. When a fully elastic and friction-free collision behavior is assumed, the particle velocities of neighboring spherical particles are only positively correlated in the Stokes regime and this correlation increases with decreasing particle-fluid density ratio (compare Fig. 5a and 5b). Inelastic collisions with additional friction lead to a significant correlation of particle velocities at higher Reynolds numbers (Re = 20 and Re = 1000) but do not change the correlation in the Stokes regime (compare Fig. 5b and 5c). This flow regime dependent influence of density ratio and collision parameters was already found for the density distribution functions and the drag forces (see Tab. 3 and Fig. 3). The corresponding results obtained for the spherocylinders (Fig. 5d) also show strong particle velocity correlations in the Stokes regime, which is even more pronounced compared to spheres. Especially, for the solids volume fraction j = 0.4 remarkably far-ranging correlations are obtained. This phenomenon can be explained by the entanglement of the cylindrical particles, which forces them to move with the same velocity in a dense particle system. The strong far-ranging correlation of particle velocities may also explain that the reduction of interphase drag is stronger for the spherocylinders than for the spheres when comparing the static to the moving case (see Sect. 3.1). Another phenomenon in Fig. 5 that should be addressed is the occurrence of slightly negative particle velocity correlations at larger particle separation distances in the cases where velocity correlations are strongly positive at small separation distances. The simulation system is designed by construction to keep the net momentum zero. Once a larger particle cluster starts to move collectively in a distinct direction with the surrounding flow, this will therefore excite remote particles to move in the opposite direction with the surrounding flow and result in a negative particle velocity correlation. Finally, the small but noticeable positive velocity correlation at small separation distances in Fig. 5a is attributed for the case Re = 1000 with j = 0.2 to the fact that particles have a high mobility (high granular temperature due to high Reynolds number), which makes it easier to reach other particles and form a small cluster that can persist for a certain amount of time due to high weight of the particles (r p /r f = 1000) and the fact that the formation of isolated clusters is facilitated by the low solids volume fraction. This can lead to the observed correlated particle motion, which is however not strongly pronounced.

Conclusions
Particle-resolved direct numerical simulations help to gain a deeper understanding for complex flow characteristics in particle-fluid flows. While PR-DNS are commonly used to study the flow through static arrangements of particles, the question of the transferability to dynamic gas-solid flows occurring in a wide range of applications is of high importance. In the present study, it is demonstrated that this transferability is not always given. Correlations for particlefluid forces obtained from the consideration of flow through static particle configurations should therefore be applied in dynamic unresolved particle-fluid simulations with caution.
The key findings of this study are as follows: -The mobility of fully elastic and friction-free particles has a significant influence on interphase drag in the Stokes regime, but not at elevated Reynolds numbers. The strong reduction of interphase drag in the Stokes regime obtained for mobile particles compared to static particles is stronger for cylindrical particles than for spherical particles, indicating a strong particle shape dependency. A lower particle-fluid density ratio further decreases the interphase drag, when considering spheres.   Fig. 1).
-Outside the Stokes regime, the collision parameters that model elasticity and friction play an important role. Inelastic collisions with friction lead to an increase of interphase drag when considering spheres. -The Reynolds number dependent influence of particle density ratio and collision parameters can also be found in the tendency for particle clustering and the correlation between the velocities of neighboring particles. -The magnitude of drift forces relative to the magnitude of drag forces is very high for both of the analyzed particle shapes, especially in dynamic systems. A static lift force correlation based on the particle orientation with respect to the slip velocity cannot predict the drift force sufficiently accurate. Since we have focused on time-and ensemble-averaged forces, structures and correlations in the present work, a following PR-DNS study has to shed light on the individual time-dependent particle forces and their relationship with the position and velocity of neighboring particles. Particularly, the origin of drift forces has to be studied more closely in order to derive rigorous drift force models. Moreover, the influence of collision parameters and particle shape was only shown exemplarily here and deserves a detailed investigation by considering a wider range of collision parameters and particle shapes.
Generally, the particle-resolved simulations can be assumed to be more accurate than particle-unresolved simulations, since they do not rely on constitutive models for the particle-fluid forces but inherently allow a direct calculation of the forces and thereby a derivation of constitutive models. This is also a major advantage over state-of-the-art experimental methods for particle tracking [32][33][34][35][36] that have been recently developed to track positions, velocities and orientations of non-spherical particles and result in macroscopic distributions of these quantities, but do not allow a measurement of individual particle-fluid forces. Due to the high computational demand of resolved simulations, unresolved simulations will nevertheless still be the only option for large-scale models in the foreseeable future. Thus, future particle-resolved DEM-CFD analyses should be accompanied by particle-unresolved DEM-CFD simulations to detect discrepancies between the two modelling www.cit-journal.com 2020   . Ensemble-and time-averaged correlation between the particle velocities of neighboring particles as a function of particle center distances under different flow conditions: a) moving case I -high particle-fluid density ratio with spheres, b) moving case II -low particle-fluid density ratio with spheres, c) moving case III -low particle-fluid density ratio with additional particle elasticity and friction with spheres, d) moving case II -low particle-fluid density ratio with spherocylinders.
stages. The development of accurate constitutive models for particle-unresolved simulation will therefore also remain an important research area that should -according to the findings of the present study -pay special attention to mobility effects and particle shape effects. Finally, as heat transfer plays a decisive role in various applications, it should also be examined how mobility effects and particle shape effects influence the heat transfer in thermal particle-laden flows and the validity of constitutive heat transfer models.