Response Behavior of Nonspherical Particles in Homogeneous Isotropic Turbulent Flows

In this study, the responsiveness of nonspherical particles, specifically ellipsoids and cylinders, in homogeneous and isotropic turbulence is investigated through kinematic simulations of the fluid velocity field. Particle tracking in such flow field includes not only the translational and rotational components but also the orientation through the Euler angles and parameters. Correlations for the flow coefficients, forces and torques, of the nonspherical particles in the range of intermediate Reynolds number are obtained from the literature. The Lagrangian time autocorrelation function, the translational and rotational particle response, and preferential orientation of the nonspherical particles in the turbulent flow are studied as function of their shape and inertia. As a result, particle autocorrelation functions, translational and rotational, decrease with aspect ratio, and particle linear root mean square velocity increases with aspect ratio, while rotational root mean square velocity first increases, reaches a maximum around aspect ratio 2, and then decreases again. Finally, cylinders do not present any preferential orientation in homogeneous isotropic turbulence, but ellipsoids do, resulting in preferred orientations that maximize the cross section exposed to the flow.


Introduction
Nowadays, the use of numerical simulation techniques to assist the development and optimization of industrial processes dealing with turbulent multiphase flow has been included as one more step in their layout. Examples of them include pneumatic conveying, fluidized bed reactors, cyclones, classifiers, or flow mixers. Industrial sectors where such processes are important are the chemical, food, or paper industries as well as electric energy production. Due to the complexity of the involved flow, a great majority of simulations are carried out under Reynoldsaveraged Navier-Stokes (RANS) in connection with an appropriate turbulence model to describe the turbulent dynamics of the carrier phase.
Two main frames are employed for the description of complex multiphase flows: the two-fluid model or Euler-Euler and the discrete particle models or Euler-Lagrange. In both of them, particles are approximated as point masses being transported in the carrier phase flow field; the solution of the flow around individual elements is usually too expensive and cannot be afforded. In the two-fluid model, both phases are conceived as two interpenetrating continua [1] whose properties are described by sets of partial differential equations. In the Euler-Lagrange approach, the discrete elements are considered as individual objects whose dynamics is governed by a Lagrangian motion equation. Therefore, to obtain the discrete phase variables in the computational domain, a large enough number of discrete element trajectories must be computed. On each particle, appropriate forces act reflecting the various microprocesses taking place at the element scale such as fluid-particle turbulent interaction, particle-(rough) wall interactions, and interparticle collisions [2]. Such technique is especially appropriate for the description of disperse multiphase flow, where usually particles have a size distribution, in confined domains where particle-wall collisions play a predominant role as pneumatic conveying, separation, and classification processes. Both techniques, twofluid model and Euler-Lagrange, have been applied mainly considering spherical particles. This means that the forces due to the flow (drag and lift) as well as the microprocesses modeling, wall-particle and inter-particle interactions, are assumed to be for spherical-shaped elements [3]. In practical situations, however, nonspherical particles are encountered, either of irregular shape, either with welldefined shapes (fibers or granulates). For example, the paper industry uses large amounts of turbulent liquid to handle and transport the fibers that compose the paper pulp. Besides, such particles in the flow experience particle Reynolds numbers larger than one, Re > 1. For such particles, the most relevant transport mechanisms such as aerodynamic transport, wall-particle interactions, and interparticle collisions are substantially different than those for spherical particles.
With the objective of performing the numerical simulation of turbulent flows laden with nonspherical particles, additional information about the forces and torques due to flow (drag and lift forces and pitching and rotational torques due to the shear flow and particle rotation) is needed. It is known that for regular nonspherical particles, that is, ellipsoids or fibers, such forces depend on particle orientation with respect to the flow. For instance, fiber orientation plays a major role in chemical processes as injection, compression molding, or extrusion in which the mechanical properties of the suspensions are determined by the orientation distribution.
For the Stokes regime, particle Reynolds number much lower than 1, the behavior of the nonspherical particles can be determined by analytical methods. Forces and torques acting on an ellipsoidal particle were analytically computed by Jeffery [4]. In a series of papers, Brenner determined the forces due to the flow acting on arbitrary-shaped nonspherical particles in the Stokes regime under different flow configuration by means of theoretical methods [5]. In the creeping flow regime, also with particle Reynolds number much lower than 1, Bläser [6] computed the forces acting of the surface on an ellipsoid in free motion for different flow situations, which allow him to suggest a simple criterion for particle breakup. The drag coefficient for particle Reynolds number higher than 1 must be obtained by experiments, physical or numerical, as the analytical methods are not applicable any more.
The experimental studies to determine the drag coefficients for nonspherical particles employ wind tunnels or sedimentation vessels. For a moderately wide particle Reynolds number range, there exist results for thin discs [7], isometric irregular particles [8], cylinders and plates [9], discs [10], and discs and cylinders [11]. Drag coefficients were developed in all cases only for certain particle orientations. Compiling such results, different correlations have been developed in terms of particle shape [12][13][14]. As representative parameter of the particle shape, two options appear to be dominant: the spherical particle equivalent diameter d p , and the sphericity, defined as the ratio between the surface of the spherical particle equivalent diameter and the actual surface of the nonspherical particle. None of such correlations takes into account the dependence of the drag coefficient with the particle orientation in the flow. There exist some correlations that consider such dependence [15][16][17] and, therefore, they are appropriated to be implemented in a Lagrangian computation scheme. Nevertheless, corresponding results for lift and pitching torque coefficients are still not equally available. One of the first results for the different coefficients in terms of the orientation of elliptic particles was obtained by Hölzer and Sommerfeld [18] using the lattice Boltzmann method (LBM). Vakil and Green [19] used DNS to study the drag and lift coefficients of cylindrical particles depending on their orientation and aspect ratio, for Reynolds number up to 40, providing a correlation for them. In such work, the underlying flow field was assumed to be uniform and the flow around the particle was completely resolved. More recently, Zastawny et al. [20] applied DNS in the frame of the implicit mirroring immersed boundary (MIB) method to obtain the flow coefficients for four different ellipsoids. The authors provide specific correlations for the drag, lift, pitching torque, and rotational torque coefficients depending on the particle Reynolds number and orientation but without including the effect of the aspect ratio. The covered Reynolds number range was up to 300. Ouchene et al. [21] determined with DNS the drag, lift, and pitching torque coefficients for prolate ellipsoids with aspect ratio up to 32 and adjusted their results to proper correlations that include the effects of particle orientation and aspect ratio up to a Reynolds number of 300.
The first numerical computations of very small nonspherical particles in pseudoturbulent flow were performed by Fan and Ahmadi [22] and Olson [23]. The hydrodynamic forces and torques were computed by the theoretical coefficients of the Stokes regime. Olson [23] estimated the time step for the translation and rotation motions in function of the fiber length, obtaining the corresponding dispersion coefficients. Fan and Ahmadi [22] showed that the dispersion of both, translation and rotation, was reduced with the fiber length. However, Olson [23] found a different result in the case of ellipsoidal particles. Lin et al. [24] investigated numerically the distribution of the orientation of the fibers in a developing mixing layer, comparing the obtained results with experiments. The fiber length was smaller than the Kolmogorov scale, so they employed the forces due to the flow of the Stokes regime. Zhang et al. [25], Mortensen et al. [26], and Marchioli et al. [27] studied the transport and deposition of ellipsoidal particles in a turbulent channel flow using direct numerical simulation (DNS). Again the hydrodynamic forces and torques were computed with Stokes regime expressions. Beyond the Stokes regime, van Wachem et al.
[28] and Ouchene et al. [29] studied a turbulent channel flow laden with ellipsoidal particles using LES and DNS, respectively, employing the flow coefficients developed by themselves in previous works.
Rosendahl group developed a model for the numerical computation of cylindrical and superellipsoidal particles in laminar and turbulent flows in the intermediate Reynolds numbers regime [30][31][32]. Particle angular velocity and orientation were computed by means of the Euler parameters. Using a linear relationship between the drag coefficient and the ellipsoid parameters, it was possible to establish a correlation valid up to Reynolds numbers of 1000. To estimate the influence of the orientation, a correlation between the maximum (90°) and minimum (0°) drag was employed. Drag force was calculated using the projected area perpendicularly to the flow. The lift force was expressed in function of the drag coefficient and particle orientation. Other lift forces, such as those due to the fluid velocity gradients or particle rotation, were not considered. The study case was a combustion chamber with straw particles, which were quite well approximated by cylinders. It was found that straw particles were better dispersed than spheres [30], a fact that properly illustrates the importance of the correct modeling of nonspherical particles motion. In a later work [32], other forces such as added mass and pressure force were also included. Drag coefficient was computed using the Ganser [15] correlation, making it possible to numerically compute the biomass combustion chamber.
This contribution aims to study the motion of nonspherical particles immersed in homogeneous isotropic turbulent (HIT) velocity fields built from kinematic simulation at moderate Reynolds numbers. Computations were performed in a tailored in-house code. Properties analyzed include the Lagrangian time autocorrelation function, the translational and rotational particle response, and preferential orientation of the nonspherical particles in the turbulent flow, all of them in terms of particle aspect ratio and inertia.

Governing equations 2.1 Coordinate systems
To build the trajectory of a regular nonspherical particle, it is necessary to solve for its translational as well as rotational motion. However, whereas translation is solved in an inertial frame, rotation is solved referred to the so-called particle frame. Thus, the relevant coordinate frames and the transformations between them have to be introduced. Figure 1 illustrates, in the case of a cylindrical particle, the employed coordinate systems: is the particle frame, whose origin is in the particle center of mass and its axes are the particle principal axes; and x 00 ¼ x 00 y 00 z 00 Â Ã is the comoving frame, which has its origin at the same point than particle frame but its axes are parallel to the inertial frame axes. In the particle frame, the z 0 axis coincides with the particle symmetry axis and its position with respect to the comoving frame determines particle orientation.
Goldstein [33] gives the transformation between the comoving and particle coordinate systems, which is frequently employed in regular nonspherical particle tracking [31].
A is the orthogonal matrix that performs the transformation. Its components are the direction cosines of the particle axes in the comoving frame, written in function of Euler angles θ; ϕ; ψ ð Þ. Such Euler angles are defined according to the x-convention of [33]: A ¼ cos ψ cos ϕ À cos θ sin ϕ sin ψ cos ψ sin ϕ À cos θ cos ϕ sin ψ sin ψ sin θ À sin ψ cos ϕ À cos θ sin ϕ cos ψ À sin ψ sin ϕ þ cos θ cos ϕ cos ψ cos ψ sin θ sin θ sin ϕ À sin θ cos ϕ cos θ (2) The time evolution of such Euler angles depends on the particle angular velocity regarding the particle frame axes. However, there is a difficulty in the sense that such time evolution equations present an unavoidable singularity. Therefore, instead of the Euler angles, the Euler parameters ε 1 ; ε 2 ; ε 3 ; η ð Þare used instead: And the transformation matrix A is written as [33]: In the present study, the initial particle orientations are assigned by means of the Euler angles. From them, the corresponding Euler parameters are computed by Eq. (3), and with them, the initial transformation matrix is evaluated using Eq. (4). The Euler parameters evolve in time following Eq. (5), where the particle angular velocities are expressed in the particle frame of reference

Particle motion equations
The nonspherical particle motion equations in a general fluid flow [34] are written as: Translational motion: Rotational motion: Here, m p is the mass of the particle, u p ¼ u px u py u pz Â Ã is the translational velocity of the particle center of mass, referred to the inertial frame, and F ¼ F x F y F z Â Ã is the external forces acting on the particle. The moments of inertia with respect to the particle frame axes are I x 0 I y 0 I z 0 Â Ã , and T x 0 T y 0 T z 0 Â Ã are the torques experienced by the particle. It should be remarked that the equations for the translation motion are computed in the inertial frame but those of the rotation motion are expressed in the particle frame. In case of the torque experienced by the particle, it has two contributions: the pitching torque, due to the noncoincidence of the particle center of mass and center of pressure (same fact that happens in an airfoil), and the rotational torque, due to the viscous resistance experienced by a rotating body inside a fluid, generated by the differences between fluid and particle rotational velocities.
In addition to a sphere, the four ellipsoids of Zawstawny et al. [20] have been chosen. They have different sphericities and aspect ratio (see Table 1). In Table 1, a denotes the major semiaxis and b the minor semiaxis.
Using DNS for ellipsoidal particles immersed in a uniform flow, Zastawny et al. [20] determined correlations for the flow coefficients (drag C D , lift C L , pitching torque C T , and rotational torque C R ). Such coefficients are written as [20]: Here, d p is the volume equivalent particle diameter or the diameter of a sphere with the same volume as the considered particle. The relative fluid velocity with respect to the particle is e u ¼ u À u p and Ω ¼ 1 2 ∇ Â u À ω p is the fluid relative rotation with ω p being the particle angular velocity.
The developed correlations depend not only on particle Reynolds number Re ¼ ρd p e u=μ and particle rotation number Re R ¼ ρd 2 p Ω j j=μ, but also on orientation φ. They are written as [20]:
Ellipsoids evaluated by Zastawny et al. [20]. a and b are the major and minor semiaxis, respectively.
Re a 8 (9) a values are listed in [20] as also coefficients b, c, r: Moreover, for the cylinders, Vakil and Green [19] developed correlations for drag and lift coefficients depending on orientation, Reynolds number based on its diameter Re D , and aspect ratio AR. In this case, such coefficients are expressed in terms of cylinder length L and diameter D: The correlations are expressed as: Coefficients κ and those γ in functions A i , i ¼ 0, 1, 2, 4: can be found in Vakil and Green [19]. However, expressions for the pitching and rotational torque coefficients are not provided in [19]. Therefore, for the cylinders, the approach of [31] has been assumed. In [31], the distance between the center of mass and the center or pressure in a cylinder, l CP , in terms of AR and φ, was proposed to be: Then, the pitching torque T P is just the cross-product between the particle orientation unitary vector and the resultant force acting on it times l CP . Nevertheless, this torque is computed in the inertial frame of reference, so it should be transformed to the particle frame before being included in Eq. (7) to calculate the particle angular velocity. The approach to compute the viscous rotational torque T R is to integrate along the particle length the torque due to the drag force with respect the particle center of mass and it is described in [31]. This torque is given directly in the particle frame.
Particle motion equations and correlations for cylinders and ellipsoids presented in this section have been implemented in an in-house code. The numerical integration of the ordinary differential equations that govern the motion of nonspherical particles has been performed by a fourth-order Runge-Kutta method, with small enough time steps to avoid numerical instabilities [35,36]. The fluid velocity field in which particles are immersed has been built by the kinematic simulation technique described in the next section. It is known that Runge-Kutta methods do not satisfy the time-reversal property, a fact that makes such methods inappropriate for integrating energy-conserving systems, for instance. However, particle equations are dissipative systems (as they include viscous drag forces) and, for them, Runge-Kutta algorithms can be used [37] provided that the time step is small enough to keep the errors bounded.

Kinematic simulation
There exist different options to calculate the Lagrangian properties in a turbulent flow. The starting point is the trajectory equation in which the position x x 0 ; t ð Þ of a particle released at point x 0 at time t = 0 is calculated solving: Here, u x; t ð Þ is the Eulerian velocity field. If it is known, it is possible to solve Eq. (18); however, finding u x; t ð Þ is not an easy task. One possibility is to work with Lagrangian statistics but then it would be needed to close the relevant Lagrangian correlations. Another option to solve Eq. (18) is to use DNS to obtain u x; t ð Þ; however, this is computationally very expensive. A much more economical alternative is the use of kinematic simulation (KS) to compute the Lagrangian characteristics of turbulent flow fields. In this technique, stochastic fluid velocity fields are constructed in such a way that their statistical properties are in agreement with those extracted from experiments or reliable DNS. The main advantage of KS is that it employs an explicit continuous formula for computing u x; t ð Þ, so it is not needed to perform interpolation of the fluid velocity field. Moreover, KS results of two particle statistics in HIT have been validated versus DNS showing good agreement [38].
The three-dimensional Eulerian velocity field to be employed in Eq. (18) is built as a series of random Fourier modes. The velocity field is solenoidal at each realization by construction. Moreover, the energy spectrum of the Fourier modes is prescribed, for example, by a power law, so the effects of small flow scales on Lagrangian statistics are directly included. Such KS velocity field is written as [39]: k n represents the n-th wave number; coefficients A n , B n are random, uncorrelated vectors perpendicular to k n , whose amplitudes are chosen according to the prescribed energy spectrum E(k) [39]. Here, the energy spectrum has been the Kolmogorov decay law of À5/3. ω n is the n-th frequency, which determines the unsteadiness of the corresponding mode; it is written proportional to the eddy-turnover time of the n-th mode: Here, λ is a parameter of order 1 that governs the unsteadiness of the velocity field. In three-dimensional HIT flows, it has been demonstrated [38] that the statistical characteristics of two-particle diffusion are independent of λ. In particular, in this work, two values of the unsteadiness parameter of 0 and 0.5 have been tested, without significant differences in the computed statistical properties. Therefore, following the suggestions of [39], the value 0.5 has been adopted in the present computations.
Because of the construction of the velocity field given by Eq. (19), k n •A n ¼ k n •B n ¼ 0, it is solenoidal trajectory by trajectory. Moreover, as shown in [40], such field includes in each realization turbulent-like patterns as eddying, straining, and streaming regions.
To validate the spherical particles tracking in the KS velocity field, the values of particle Reynolds stresses (RS) in HIT have been selected. In this configuration, Hyland et al. [41] demonstrated that, as the fluid turbulence is homogeneous, particle RS only depend on time and they can be written as u 0 ð Þδ ij , that is, they are an isotropic tensor. Moreover, in the asymptotic limit, q t ! ∞ ð Þ¼βT L = 1 þ βT L ð Þ , where β is the inverse of particle relaxation time (see Eq. (21) below) and T L is the Lagrangian time scale of fluid turbulence. Figure 2 presents the numerical results for particle RS computed with KS and the asymptotic expression q t ! ∞ ð Þ. As Figure 2 readily shows, the asymptotic particle RS are very well reproduced by the numerical particle tracking in the KS velocity field in the range of two decades for βT L .

Numerical simulation
Computations were performed in a tailored in-house code. The turbulent velocity field generated with KS resembles one of the fields worked in [39]. Such velocity field is characterized by a fluctuating velocity u' = 1 m/s, a fluid Reynolds number of 10 4 resulting in a Kolmogorov length scale η K ≈ 6:286 mm, associated Kolmogorov time scale τ K ≈ 10 ms, and a fluid integral Lagrangian time scale of turbulence T L ¼ 0:56 s. Those values are matched by the present KS.
The regular nonspherical particles studied have been the ellipsoids in [20] and the cylinders in [19]. In all cases, particles have the same particle volume equivalent diameter d p ¼ 200 μm, hence much smaller than η K . Therefore, such particles can be thought as immersed in a uniform flow field. The Stokes number has been modified by adjusting the material density of particles being the Stokesian particle relaxation time defined as: If the Kolmogorov time scale τ K is taken as the fluid time scale, particle Stokes number is defined as St ¼ τ p =τ K . According to this nondimensional number, three particle inertia classes are considered: light (St ≈ 0:5), intermediate (St ≈ 10), and heavy (St ≈ 100). However, as cases with Re > 1 are considered, an effective particle relaxation time is introduced as τ p, eff ¼ τ p =ReC D , allowing the introduction of an effective Stokes number St eff ¼ τ p, eff =τ K . Therefore, the values of such effective Stokes number are St eff ≈ 0.3 (light), 5 (intermediate), and 40 (heavy).
Simulations proceed in the following way: for each KS realization of HIT fluid velocity field, a particle is located in the center of the domain with zero initial velocity; particle translational and rotational motion is computed from Eqs. (6) and (7), its orientation is calculated from Eq. (1), and its trajectory is built; particle tracking lasts for around 10 fluid integral time scales; and particle properties are stored every second for evaluation. Such process is carried out a sufficient number of times to reach significant statistical results. In this study, statistics has been performed based on 10 5 KS realizations.
In the following section, the results of the particle Lagrangian time autocorrelation function, the translational and rotational particle response, and preferential orientation of the nonspherical particles in the turbulent flow are analyzed as function of their shape and effective Stokes number.

Results and discussion
The Lagrangian autocorrelation function R L, t τ ð Þ for translational motion is expressed as: τ represents the time delay. With the objective of making results independent of particle injection conditions, statistics are started to be collected after 2 s. The obtained results for the Lagrangian autocorrelation function (LAF) of the ellipsoids of Zastawny et al. [20] are presented in the left part of Figure 3, including the results for spherical particles, whereas those of the cylinders of Vakil and Green [19] are in the right side of such figure. Horizontal axis is the nondimensional time delay, τ=T L . As the lighter particles have an autocorrelation function nearly equal to that of the fluid (tracer limit), the corresponding curves are not shown in Figure 3. Therefore, only the curves for intermediate and high inertia particles are presented. Moreover, in Figure 3 also the fluid Lagrangian (brown curve) and Eulerian (cyan curve) R L, t τ ð Þ's are included for comparison. As it can be readily seen from Figure 3, higher inertia particles are characterized by larger integral Lagrangian time scales (ILTSs) (defined as the integral up to infinity of R L, t τ ð Þ), as a result of their smaller responsiveness to the turbulent fluctuations. Same as in [23], all particle LAFs are mainly in between the fluid LAF and Eulerian autocorrelation function (EAF). As a consequence of inertia, for the smallest values of τ, the heavy R L, t τ ð Þ overcomes the fluid EAF, differently from [23] who considered only noninertial particles.
Moreover, for intermediate inertia, the curves for all particle shapes nearly collapse in a single curve. On the other hand, a shape effect is noticeable for the heaviest particles, where the various shapes present differences in their curves. It is interesting to realize that ILTSs of higher aspect ratio (AR) are below those of smaller AR, for both ellipsoids and cylinders. This effect is a Reynolds number effect due to the dependence of drag coefficient on shape and AR: an interaction between translation and rotation motions occurs that results in a spreading of the particle effective Stokes number. As a consequence, particles with higher Reynolds numbers also have larger effective inertia (reflected on an increased Stokes number) and, therefore, their LAF decreases slower, implying a higher ILTS. In the ellipsoids case, it happens that those of typ. 2 present a R L, t τ ð Þ curve slightly over that of the spherical particle, as they have lower effective Stokes number. Also, the LAF curve for the disc-like particles in this case is very similar to that of the fiber.
In an analogous way to translational LAF, a rotational autocorrelation function (RAF) R L, r τ ð Þ can be defined in terms of the time delay τ as: where the particle angular velocity is denoted by ω p . The obtained results for the Lagrangian rotational autocorrelation function of the ellipsoids of Zastawny et al. [20] are shown in the left part of Figure 4, whereas those of the cylinders of Vakil and Green [19] are in the right side of such figure. Again, horizontal axis is the nondimensional time delay, τ=T L . Similar to the case of translational motion, the angular velocities of heavy particles keep correlated for longer times than those of lighter particles. Such correlation time for ellipsoidal particles is much shorter than that of the translational motion. Also, for all inertia cases, the RAF of disc-like particles drops quicker than for the prolate ellipsoids. As mentioned for the translational correlations, the RAF curves for the prolate ellipsoids collapse for the lighter particles, but they show noticeable differences for the intermediate and large inertia particles demonstrating an effect of the aspect ratio on R L, r τ ð Þ. Ellipsoid 2, with the smaller AR, has the higher RAF curve of all prolate ellipsoids, while Ellipsoid 1 and the fiber have very similar rotational correlation functions.
In the case of cylinders (Figure 4, right), the RAF curves for all AR and inertias are different. For the smallest inertia particles, RAF decreases with increasing AR, similar to what was found for LAFs in the translational motion. Moreover, the R L, r τ ð Þ curve presents negative values for the two largest aspect ratios of 10 and 20. For the intermediate particles, the RAF curves keep the same decreasing trend with increasing aspect ratio as the light particles; however, in this case, correlation times for angular velocities are very much increased and they are significantly higher than for ellipsoids. The previous trend is reversed for the heavy particles as the RAF curves augment with increasing AR; however, as Figure 4 (right) suggests, an asymptotic value for L/D seems to exist because the curves for AR = 10 and 20 are very close to each other. The change of behavior of the RAF with increasing inertia could be due to the fact that for the small and intermediate inertia, the particle relaxation time for rotation reduces with growing AR, whereas for heavy particles, such relaxation time behaves in the opposite way. Nevertheless, this fact must be further investigated, possibly using fully resolved simulations.
Next, the response of the nonspherical particles to the fluid fluctuating velocities is analyzed for both translation and rotation motions. Figure 5(a) shows the behavior of the particle's relative linear root mean square (rms) velocity, that is, u 0 p =u 0 , where u 0 is the fluid rms fluctuating translational velocity and u 0 p denotes the same quantity but for the particles. The aspect ratio is in the horizontal axis, whereas the different curves correspond to the various inertia cases. On the one side, as it could be anticipated, u 0 p reduces with increasing particle inertia because the more inertial particles are not able to follow all fluid velocity fluctuations. In fact, as it was found for the LAF, the less inertial particles present, for both ellipsoids and cylinders and for all values of AR, the same fluctuating velocities as the fluid, indicating that they behave as fluid tracers.
As inertia increases, u 0 p =u 0 decreases monotonically, as expected. However, it increases with growing aspect ratio for both cylinders and ellipsoids. There is one exception that spherical particles have values of u 0 p =u 0 slightly above those of the Computed R L, r τ ð Þ curves for ellipsoidal particles [20] (left) and cylindrical particles [19] (right). ellipsoid with AR = 1.25. The trend of increasing particle fluctuating velocities with aspect ratio is consistent with the aforementioned fact that effective Stokes number tends to reduce with growing AR; therefore, particles with lower AR respond less to the fluid fluctuations than the more elongated ones. This result has been obtained for particles much smaller than Kolmogorov length scale; therefore, such particles can be considered to be immersed in a uniform flow field, just in the same conditions as the flow coefficient correlations were developed. In the study of Hölzer and Sommerfeld [42], nevertheless, a different result was found. Using a DNS based on the lattice Boltzmann method (LBM), Hölzer and Sommerfeld [42] obtained that relative particle fluctuating velocity reduced with increasing AR. The main difference with the present work is that Hölzer and Sommerfeld [42] employed particles with size well above η K . The authors explained the fact arguing that particles averaged the fluid fluctuations on their surface, which augmented with increasing AR. Let us remark that both results are not conflicting as the size range of the employed particles in the two studies is very different. Further work combining DNS with nonspherical point particles smaller than η K is necessary to explain this point. Figure 5(b) presents the behavior of the particle angular rms velocity, ω 0 p . In homogeneous and isotropic turbulence, the angular velocity of spherical particles is zero because of viscous damping and the absence of pitching torque. The situation is different in nonspherical particles because in them the geometrical and pressure centers do not coincide, so a net pitching torque is produced that promotes nonzero angular velocities. As illustrated in Figure 5(b), ω 0 p increases with AR and reaches up to a maximum of AR ≈ 2 and decreases with higher values of aspect ratio. The shape of the curve is the same for ellipsoids and cylinders. Such behavior was also found in [42], and it was explained observing that, with increasing AR, the moment of inertia along the major axis reduces and along the minor axis increases, which would lead to higher and lower ω 0 p , respectively. On the other hand, and similar to what happened with u 0 p , for inertial particles, ω 0 p decreases as inertia augments. In the following, the correlation relative velocity direction-particle orientation is analyzed depending on inertia and aspect ratio. A well-known fact is that regular nonspherical particles falling through a still liquid at intermediate Reynolds numbers tend to be oriented in a determined direction. Cylinders and prolate ellipsoids are prone to keep their symmetry axis (z' in Figure 1) perpendicular to the flow, thus maximizing drag. Differently, discs and oblate ellipsoids tend to move with the symmetry axis aligned with the flow, also maximizing drag [43]. However, spheroidal Stokes particles only show a preferential orientation if a persistent velocity gradient exists [27]. Therefore, in HIT flow where there are no mean velocity gradients, a Stokes particle will not have any preferred orientation.
Newsom and Bruce [44] analyzed the influence of turbulence on preferential alignment of quite elongated fibers with Re ≈ 1. As explained by [44], preferential orientation of such fibers falling through a still fluid can only be clarified if fluid inertial effects are considered. In the Stokes regime, Khayat and Cox [43] demonstrate that the force distribution on the fiber is symmetrically distributed along its axis, independent of its orientation regarding the flow and, as a result, the fiber experiences a zero net torque. Beyond the Stokes regime, Re > 1, when the fiber has an oblique orientation with respect the flow, such distribution of the force is not any more symmetric and it experiences a net pitching torque. Such torque will promote a rotation that drives the fiber to adopt an orientation where its symmetry axis is orthogonal to the relative flow. Interestingly, if the fiber is oriented orthogonal or parallel to the flow, the net experienced torque is zero, due to the symmetry of the force distribution; however, in the first case, this situation is stable, while in the second case, where the centers of gravity and pressure do not coincide, this situation is unstable.
Previous reasoning is valid too for another kind of nonspherical shapes as disclike, cylindrical, or ellipsoidal [45]. For high Reynolds numbers, that is, Re > 100 appears a secondary motion overimposed to the particles predominant movement direction. Such secondary motion is promoted by a wake instability and vortex detachment from the rear surface of the particles. Two main kinds of secondary motion can be observed: large quasi-periodic swings along the main path, and a more or less chaotic tumbling forming a definite angle with the main motion direction. There is a coupling between this kind of oscillatory motion and the wake instability [10]: a vortex detachment follows at the end of a particle swing. Nevertheless, in the present study, such secondary motions do not appear as the considered particle Reynolds number is not large enough, that is, Re < 40.
Let now θ be the angle formed by the relative velocity, u À u p , and the particle symmetry axis, z'. Therefore, cos θ can be used to determine the orientation of the nonspherical particle with respect to the relative flow. In this work, particle preferential orientation is determined computing cos θ j jalong the trajectories of 10 5 particles. Computed values of cos θ j jare sorted in equally distributed bins between 0 (particle axis orthogonal to relative velocity) and 1 (alignment between particle axis and relative velocity), and the corresponding probability density functions (Pdfs) are determined. Such Pdfs are shown in Figure 6 in terms of cos θ j j. Figure 6(a) shows the results for the prolate ellipsoids in terms of particle inertia, Figure 6. Orientations probability density functions (Pdfs) of prolate ellipsoids (a) and cylinders (b) regarding the relative flow direction. and Figure 6(b) presents the curves for the cylinders also depending on their inertia. Each inertia class is plotted in a separated frame. For the discs case, results are presented in Figure 7.
As it is observed in Figure 6(a), it is found that prolate ellipsoids do manifest preferential orientation with respect to the relative velocity. Of course, spherical particles do not have a preferred orientation and the corresponding Pdf is a horizontal line (black color). Prolate ellipsoids have a preference for orientating its symmetry axis orthogonal to the relative flow, tending to maximize the drag, similar to what occurs in particle sedimentation studies. The orientation preference increases with inertia, which is quite similar for all aspect ratios considered in this study.
On the other hand, as it is presented in Figure 6(b), cylinders seem not to have any preferred orientation in the HIT KS velocity field, being all the curves pretty flat. Only for the case of higher AR and lowest inertia, the curves show a trend to be slightly higher for values of cos θ j jcloser to one than to zero. Such result is qualitatively similar to the DNS computations of [27] in the central region of the channel.
For the discs, Figure 7 shows that there is a clear trend of the particle symmetry axis to be aligned with the relative flow, again maximizing drag, similar to the results obtained for sedimenting particles in stagnant fluid. Such trend is more marked when particle inertia increases.

Conclusions
In this study, regular nonspherical particle responsiveness to HIT flows has been investigated in combination with KS of fluid velocity field. The main results obtained are the following: the particle LAF reduces when particle AR is augmented, because effective particle inertia decreases if aspect ratio increases; this is true for both translational and rotational time autocorrelation functions. In the case of cylinders, R L, r τ ð Þ is much higher than for ellipsoids, a fact that requires further clarification through particle resolved simulations. Additionally, the fluctuating particle velocity increases for growing AR in the considered case of particles much smaller than Kolmogorov length scale; such behavior is contrary, although not conflicting, to the findings of [42] for fully resolved particles with sizes larger than the Kolmogorov length scale. For both ellipsoids and cylinders, the particle angular rms velocity first increases with aspect ratio, reaches a maximum of AR ≈ 2, and then decreases again, which is explained because with increasing aspect ratio, the moment of inertia around the longitudinal axis decreases and around the radial axis increases, which would lead to higher and lower rms angular velocities, respectively. Finally, in agreement with Marchioli et al. [27], cylinders seem not to prefer any specific orientation in the KS HIT velocity field; however, prolate ellipsoids tend to be oriented with its symmetry axis orthogonal to the relative flow, maximizing the drag. Oblate ellipsoids and disc-like particles also show a preferential orientation, tending to align their symmetry axis with the relative flow velocity.