Stochastic dispersion of ellipsoidal fibers in various turbulent fields

Author's Accepted Manuscript

Stochastic dispersion of ellipsoidal fibers in various turbulent fields M.M. Tavakol, O. Abouali, M. Yaghoubi, G. Ahmadi

www.elsevier.com/locate/jaerosci

PII: DOI: Reference:

S0021-8502(14)00174-8 http://dx.doi.org/10.1016/j.jaerosci.2014.10.008 AS4834

To appear in:

Journal of Aerosol Science

Received date: 9 July 2014 Revised date: 29 October 2014 Accepted date: 30 October 2014 Cite this article as: M.M. Tavakol, O. Abouali, M. Yaghoubi, G. Ahmadi, Stochastic dispersion of ellipsoidal fibers in various turbulent fields, Journal of Aerosol Science, http://dx.doi.org/10.1016/j.jaerosci.2014.10.008 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting galley proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Stochastic dispersion of ellipsoidal fibers in various turbulent fields M.M. Tavakol1, O. Abouali1, , M. Yaghoubi1,2, G. Ahmadi3 1

School of Mechanical Engineering, Shiraz University, Shiraz, Iran 2

3

Academy of Sciences, I. R. IRAN

Department of Mechanical and Aeronautical Engineering, Clarkson University, Potsdam, NY, USA

Abstract In this study dispersion and deposition of fibers of various sizes in fully developed turbulent pipe and duct flows were analyzed for a range of flow Reynolds numbers. Fibers translational and rotational equations of motion were solved using the Lagrangian approach assuming one-way interaction. The influences of turbulent fluctuations were included using appropriate stochastic models for fluctuating velocity and velocity gradients components. Performance of the stochastic model was evaluated for several test cases including dispersion of fibers in a synthesized isotropic homogenous turbulent flow field, inhomogeneous turbulent flow field, and fully developed turbulent pipe and duct flows. The simulation results were compared with earlier numerical and experimental studies for flow statistics and deposition velocity and good agreements were observed. It was shown that using appropriate stochastic models leads to satisfactory evaluation of dispersion and deposition of fibers in turbulent flows.

Keywords: Turbulent flow, Fibers, Stochastic modeling, Velocity gradient, Particle deposition, Eulerian-Lagrangian method



Corresponding author: O. Abouali, School of Mechanical Engineering, Shiraz, Iran. E-mail address: [email protected], [email protected]



1

1. Introduction The physics of particle transport and deposition is important in many areas of fluid engineering. Fouling in heat exchangers, deposition of particulate pollutant in respiratory tract, transport and deposition of particles in microelectronics and paper industries are some examples. In the recent years, commercial CFD software for simulation of particle transport and deposition were developed, which typically assume that the particles are spherical. However, majority of the particles suspended in airflows are non-spherical and many are elongated fibers with large aspect ratios. The principal factors affecting transport and deposition pattern of fibers include aerodynamic properties of these particles, their size and aspect ratio, and the flow velocity field. Several studies related to particle motion were reported in the literature. Using the concept of free flight, Friedlander and Johnstone (1957) developed a theory for particle deposition in turbulent flows. Ahmadi (1970) reported analytical and numerical simulations of dispersion of small suspended particles in turbulent flows. Reeks (1977) and Maxey (1987) studied the effect of crossing trajectories and gravitational settling on dispersion of particles. Ounis and Ahmadi (1989, 1990, 1991) conducted numerical simulations and reported theoretical models for diffusivity of particles in isotropic and turbulent shear flows. They noted the important role of shear-induced lift force on particle dispersion across a shear field. Cleaver and Yates (1975) studied the particle deposition mechanism due to turbulent burst and inrush. Further progress along this line were reported by Fichman et al. (1988), and Fan and Ahmadi (1993, 1994, 1995a, 1995b, 2000). Numerical simulations of particle deposition in laminar and turbulent duct flows were performed by Li and Ahmadi (1993), He and Ahmadi (1999), Zhang et al. (2000, 2001), and Tian and Ahmadi (2007). More recently, Tian et al. (2012) had performed experimental and numerical studies to evaluate the motion of ellipsoidal fibers in low Reynolds number flows. Early models for particle deposition were based on the stopping distance concept. According to this model, when particles reach their stopping distance, they are transported to the wall by the so-called free flight mechanism. Cleaver and Yates (1975) reported large discrepancies of the free flight model predictions with the experimental data for turbulent flow regimes.



2

In the last decade, the Direct Numerical Simulation (DNS) of Navier-Stokes equations provided a powerful tool for studying particle dispersion and deposition in turbulent flows at moderate Reynolds numbers. One-way coupled DNS of fluid flow with a Lagrangian particle tracking has been used to obtain detailed information on particle wall interactions for turbulent duct flows. McLaughlin (1989), Brooke et al. (1992), Squires and Eaton (1991), Soltani and Ahmadi (1995, 2000) and Soltani et al. (1998) used DNS to study particle deposition rate in turbulent channel flows. Although, DNS provides an accurate description of particle laden flows, it is restricted to low to moderate Reynolds numbers and simple geometry ducts. Some of the recent studies for the fiber deposition in complex respiratory passages were performed for laminar flows (Dastan et al. 2014). Recently, Tian and Ahmadi (2013) reported fiber transport and deposition in the upper airways using the Reynolds stress transport model. DNS or fiber transport and deposition in turbulent duct flows were reported by Zhang et al. (2001, 2007), and Soltani et al. (2000). The presented literature survey shows that most of the earlier studies were concerned with the dispersion analysis of spherical particles. For non-spherical particles, earlier studies were reported for laminar shear flows, and there is only limited number of studies regarding turbulent dispersion of non-spherical particles.

The aim of the present study is to develop proper

stochastic models for analyzing dispersion and deposition of ellipsoidal particles in isotropichomogeneous and inhomogeneous turbulent flows. In the following sections, first the mathematical background of the fiber motion is described. This is followed by the discussion of the method to incorporate the effects of turbulence fluctuations on fiber dispersion. Finally, performance of the developed model for dispersion of fibers in various turbulent flows are evaluated and discussed.

2. Mathematical background 2.1. Creeping flow formulation for fiber motion The motion of a non-spherical particle is described using translational and rotational equations of motion. Translational motion is described using the equations of balance of linear momentum for the particles .Rotational motion is described by the equations of balance of angular momentum. To evaluate the motion of an ellipsoidal fiber, three Cartesian coordinate



3

systems should be defined as shown in Fig.1. The (x, y, z) coordinate system is the inertial coordinate system, while (x', y', z') coordinate system is the particle frame with the origin located at the fiber center of mass with the z' along the major axis of the fiber. The third coordinate system (x", y", z") is the fiber co-moving frame with its origin at the particle centroid and its axes parallel to the inertial coordinate system. The inertial coordinate system is used to describe translational equations of motion while the rotational equations of motion are considered in the particle frame. The transformation from one to other coordinate system can be performed by means of the Euler angles , , . Goldstein (1980) suggested three successive rotations by means of Euler angles to transform between different coordinate systems. Transformation between co-moving coordinate system and particle coordinate system is given as, ‫ ݔ‬ᇱ ൌ ‫ ݔܣ‬ᇱᇱ

(1)

where A denotes the transformation matrix (Goldstein,1980): …‘• ɗ …‘• ߠ െ …‘• ߠ •‹ ߮ •‹ ߰ െ ‫ ܣ‬ൌ ൥ •‹ ߰ …‘• ߮ െ …‘• ߠ •‹ ߮ …‘• ߰ •‹ ߠ •‹ ߮

…‘• ߰ •‹ ߮ ൅ …‘• ߠ …‘• ߮ •‹ ߰ െ •‹ ߰ •‹ ߮ ൅ …‘• ߠ …‘• ߮ …‘• ߰ െ •‹ ߠ …‘• ߮

•‹ ߰ •‹ ߠ …‘• ߰ •‹ ߠ൩(2) …‘• ߠ

Due to a possible singularity in the inversion of the transformation matrix the quaternions (1, 2, 3, ) are used (Goldstein, 1980; Fan and Ahmadi, 1995; Dastan and Abouali, 2013).

ʹሺߝଵ ߝଶ ൅ ߝଷ ߟሻ ʹሺߝଵ ߝଷ െ ߝଶ ߟሻ ͳ െ ʹሺߝଶ ଶ ൅ ߝଷ ଶ ሻ ଶ ଶ ‫ ܣ‬ൌ ቎ ʹሺߝଵ ߝଶ െ ߝଷ ߟሻ ͳ െ ʹሺߝଵ ൅ ߝଷ ሻ ʹሺߝଷ ߝଶ ൅ ߝଵ ߟሻ ቏ ʹሺߝଵ ߝଷ ൅ ߝଶ ߟሻ ʹሺߝଷ ߝଶ െ ߝଵ ߟሻ ͳ െ ʹሺߝଶ ଶ ൅ ߝଵ ଶ ሻ

(2)

The quaternions are related to the Euler angles. That is, ߝଵ ൌ …‘•

ఝିట ଶ

•‹

ఏ

ߝଶ ൌ •‹

߮െ߰ ߠ •‹ ʹ ʹ

ߝଷ ൌ •‹

߮൅߰ ߠ …‘• ʹ ʹ

ߟ ൌ …‘•

(3)

ଶ

߮൅߰ ߠ …‘• ʹ ʹ

The Euler four parameters are subjected to the constraint:



4

ߝଵ ଶ ൅ ߝଶ ଶ ൅ ߝଷ ଶ ൅ ߟଶ ൌ ͳ

(4)

The translational equations of motion of an ellipsoidal fiber with mass ݉௣ are solved in the inertial coordinate system and the rotational equations of motions are expressed in the particle frame: ௗ࢜

݉௣ ‫ܫ‬௫ ᇲ ‫ܫ‬௬ ᇲ ‫ܫ‬௭ᇱ

ൌ ࡲࢎǤࢊ ൅ ࡲࢍ ൅ ࡲࡸ

ௗ௧

ௗఠೣᇲ

ൌ σ௜ ܶோ௫ ᇲ ǡ௜ ൅ ߱௬ᇲ ߱௭ ᇲ ൫‫ܫ‬௬ᇲ െ ‫ܫ‬௭ ᇲ ൯

(6)

ൌ σ௜ ܶோ௬ᇲ ǡ௜ ൅ ߱௫ᇱ ߱௭ ᇲ ሺ‫ܫ‬௭ ᇲ െ ‫ܫ‬௫ ᇲ ሻ

(7)

ൌ σ௜ ܶோ௭ ᇲ ǡ௜ ൅ ߱௬ᇲ ߱௫ ᇲ ൫‫ܫ‬௫ ᇲ െ ‫ܫ‬௬ᇲ ൯

(8)

ௗ௧ ௗఠ೤ᇲ ௗ௧ ௗఠ೥ᇲ ௗ௧

(5)

Here ‫ܫ‬௫ ᇲ , ‫ܫ‬௬ᇲ and ‫ܫ‬௭ᇱ denote the fiber principal moments of inertia: ଵ

‫ܫ‬௫ ᇲ ൌ ‫ܫ‬௬ᇲ ൌ ݉௣ ܽଶ ሺͳ ൅ ߚଶ ሻ

(9)

ହ

ଶ

‫ܫ‬௭ ᇲ ൌ ݉ ௣ ܽ ଶ

(10)

ହ

For ellipsoids the semi-minor axis of ellipsoid is a and the semi major axis is b as shown in ௕

Fig.2. The aspect ratio ߚ ൌ is the ratio of semi major axis to semi minor axis of the ellipsoid. ௔

In Eqs. (6)-(8), ߱௫ᇱ ǡ ߱௬ᇲ and ߱௭ ᇲ are the fiber rotational velocities. To update the particle orientation at each time step the time rate changes of Euler parameters are used, That is, ௗ ߝͳ

‫ ۍ‬ௗ௧ ‫ې‬ ߟ߱‫ݔ‬Ԣ െ ߝ͵ ߱‫ݕ‬Ԣ ൅ ߝʹ ߱‫ݖ‬Ԣ ‫ۍ‬ ‫ې‬ ‫ێ‬ௗߝʹ‫ۑ‬ ‫ ێ‬ௗ௧ ‫ ۑ‬ଵ ‫ݔ߱ ͵ߝ ێ‬Ԣ ൅ ߟ߱‫ݕ‬Ԣ െ ߝͳ ߱‫ݖ‬Ԣ ‫ۑ‬ ‫ێ‬ௗߝ͵‫ ۑ‬ൌ ଶ ‫ ێ‬െߝʹ ߱ Ԣ ൅ ߝͳ ߱ Ԣ ൅ ߟ߱ Ԣ ‫ۑ‬ ‫ݔ‬ ‫ݕ‬ ‫ۑ ݖ‬ ‫ ێ‬ௗ௧ ‫ۑ‬ ‫ێ‬ ‫ ێ‬ௗߟ ‫ۑ‬ ‫ۏ‬െߝͳ ߱‫ݔ‬Ԣ െ ߝʹ ߱‫ݕ‬Ԣ െ ߝ͵ ߱‫ݖ‬Ԣ ‫ے‬ ‫ ۏ‬ௗ௧ ‫ے‬

(11)

The location of particle mass center is evaluated at each time step by integrating the velocity. That is, ௗ࢞ ௗ௧

ൌ ࢜

(12)

where࢞and ࢜are location and velocity vectors in the inertial coordinate system. Assuming that fibers are ellipsoids of revolution, the corresponding forces and torques are evaluated. The hydrodynamic drag force Fh.d, for an ellipsoidal particle was derived by Brenner (1963): 

5

ࡲࢎǤࢊ ൌ ߨߤܽࡷᇱᇱ ൉ ሺ࢛ െ ࢜ሻ

(13)

where࢛is the fluid velocity vector at the fiber centroid, ߤis the fluid viscosity and ࡷᇱᇱ is the fiber resistance tensor which is defined as: ࡷᇱᇱ ൌ  ࡭ି૚ ൉ ࡷᇱ ൉ ࡭

(14)

In this equation A is transformation matrix, given by Eq. (3), and ‫ ܭ‬ᇱ is a diagonal matrix with its element given as: ࡷᇱ ࢞ᇲ ࢞ᇲ ൌ ࡷᇱ ࢟ᇲ ࢟ᇲ ൌ ͳ͸ሺߚଶ െ ͳሻȀ ൤ ࡷᇱ ࢠᇲ ࢠᇲ ൌ ͺሺߚଶ െ ͳሻȀ ൤

ଶఉమ ିଵ ඥఉమ ିଵ

ଶఉమ ିଷ ඥఉమ ିଵ

݈݊൫ߚ ൅ ඥߚଶ െ ͳ൯ ൅ ߚ൨

݈݊൫ߚ ൅ ඥߚଶ െ ͳ൯ െ ߚ൨

(15)

(16)

The gravity force is given by ࡲࢍ ൌ  ߩ௣ ܸࢍ

(17)

ସ

whereܸ ൌ ߨܽଷ ߚ is the fiber volume. ଷ

For the lift force acting on a particle in a simple shear flow, Harper and Chang (1968) presented a correlation given as, ࡲࡸ ൌ

గమ ఓ௔మ డ௨ೣ Τడ௬ భ ఔమ

భ

ȁడ௨ೣ Τడ௬ȁమ

ሺࡷᇱᇱ ൉ ࡸ ൉ ࡷᇱᇱ ሻ ൉ ሺ࢛ െ ࢜ሻ

(18)

where ͲǤͲͷͲͳ ͲǤͲ͵ʹͻ ͲǤͲͲ ‫ ܮ‬ൌ ൥ͲǤͲͳͺʹ ͲǤͲͳ͹͵ ͲǤͲͲ ൩ ͲǤͲͲ ͲǤͲͲ ͲǤͲ͵͹͵

(19)

and the lift is perpendicular to the flow direction. Other forces acting on the fibers are Brownian force, virtual mass force, pressure gradient force and Basset force. Considering the order of importance of these forces as described by Fan and Ahmadi (1995) and Kleinstreuer and Feng (2013), in the current study only the hydrodynamic drag, gravity and lift forces are retained and other forces are neglected. The hydrodynamic torques acting on an ellipsoidal fiber in the simple shear flow as derived by Jeffery (1922) is used in the present analysis. Accordingly, the hydrodynamic torques are given as,



6

ܶோǡ௫ ᇲ ൌ

ܶோǡ௬ᇲ ൌ

ܶோǡ௭ ᇲ ൌ

ଵ଺గఓ௔య ఉ ଷሺఉమ ఊబ ାఉబ ሻ

ଵ଺గఓ௔య ఉ ଷሺఉమ ఊబ ାఈబ ሻ

ሾሺͳ െ ߚଶ ሻ݀௭ ᇲ ௬ᇲ ൅ ሺߚଶ ൅ ͳሻ൫‫ݓ‬௭ ᇲ ௬ᇲ െ ߱௫ ᇲ ൯ሿ

(20)

ሾሺߚଶ െ ͳሻ݀௫ ᇲ௭ ᇲ ൅ ሺߚଶ ൅ ͳሻ൫‫ݓ‬௫ ᇲ௭ ᇲ െ ߱௬ᇲ ൯ሿ

(21)

ଷଶగఓ௔య ఉ ൣ൫‫ݓ‬௬ᇲ ௫ ᇲ ଷሺఉబ ାఈబ ሻ

െ ߱௭ ᇲ ൯൧

(22)

where the non-dimensional parameters ߙ଴ , ߚ଴ and ߛ଴ are presented by Gallily and Cohen (1979) to be: ߙ଴ ൌ ߚ଴ ൌ

ߛ଴ ൌ െ

ఉమ ఉమ ିଵ

ଶ ఉమ ିଵ

൅

െ ሺఉమ

ఉ ଶሺఉమ ିଵሻయȀమ

ఉ ିଵሻ

యȀమ ݈݊ ൤

݈݊ ൤

ఉିඥఉమ ିଵ

൨

(23)

ఉାඥఉమ ିଵ

ఉିඥఉమ ିଵ

൨

(24)

ఉାඥఉమ ିଵ

In Eqs. (20)-(22) ݀௜ ᇲ௝ ᇲ is the deformation rate tensor and ‫ݓ‬௜ ᇲ௝ ᇲ is the spin tensor which are written in the particle frame: ଵ డ௨ ᇲ

డ௨ೕᇲ

ೕᇲ

೔ᇲ

݀௜ ᇲ ௝ ᇲ ൌ ൬ ೔ ൅ ൰ ଶ డ௫ డ௫ ଵ డ௨ ᇲ

డ௨ೕᇲ

ೕᇲ

೔ᇲ

(25)

‫ݓ‬௜ ᇲ ௝ ᇲ ൌ ൬ ೔ െ ൰ ଶ డ௫ డ௫

(26)

The transformation for velocity gradient from co-moving frame to particle frame is given by: ࢛ࣔ࢏ᇲ ࣔ࢞࢐ᇲ

ൌ ࡭

࢛ࣔ࢏ᇲᇲ ࣔ࢞࢐ᇲᇲ

࢛ࣔ࢏ᇲ

where

ࣔ࢞࢐ᇲ

࡭ି૚

and

࢛ࣔ࢏ᇲᇲ ࣔ࢞࢐ᇲᇲ

(27)

are the fluid velocity gradient with respect to particle and co-moving frame.

It should be noted that when the formulation of Jeffery (1922) for hydrodynamic torques are applied, the underlying assumption is that the flow near the particle (on submicron scale) can be approximated as linear shear flow. Using Eqs. (5)-(8) it is possible to evaluate the particle position and orientation at each time step. A fourth order Runge-Kutta method is used for solving the governing equations of fiber motions as described by Eqs. (5)-(8).



7

2.2. Turbulence fluctuations modeling To be able to properly evaluate dispersion and deposition of ellipsoidal fibers in a fully developed turbulent pipe and duct flows, the effect of instantaneous turbulent velocity and velocity gradient fluctuations should be included. While there have been a number of stochastic approaches for simulating turbulence fluctuations, the simulations of instantaneous velocity gradient have been scarce. Stochastic dispersion of non-spherical particles was performed by Kvasnak et al. (2002) in connection with the formation of spray from a liquid. They examined the dispersion of deformed ellipsoidal droplets in a spray using a simulated turbulent flow field. The instantaneous fluctuating velocity and velocity gradient components at the center of each droplet were calculated with the use of probability density function (PDF) based Langevin equation. In the present study a similar procedure is used to generate instantaneous fluctuating velocity and velocity gradients, which are required for analyzing the translational and rotational motion of fibers. Most available commercial softwares use the discrete random walk (DRW) model for stochastic modeling of the turbulent velocity field. This model initially developed by Hutchinson et al. (1971) and then extended by Gosman and Ioannides (1983). In this approach the turbulent dispersion of particles is modeled as a succession of separate interactions between particles and discrete eddies. However, several investigators reported the shortcomings of such a model. Accordingly, in this approach the simulated fluctuation velocities for different direction are stochastically uncorrelated even if the non-isotropy of velocity variances can be incorporated in the generated velocity components (Mehel et al., 2010). Another option for stochastic modeling of velocity is the continuous random walk (CRW) model which is used previously by Tian and Ahmadi (2007) among others. The CRW approach is constructed based on the Langevin equation for description of Brownian motion. The governing equations for classical Langevin-equation in Cartesian coordinate system are given as,

݀‫ݑ‬ଵ ൌ െ



௨భ ்ಽ

݀‫ ݐ‬൅ ට

തതതതത మ ଶ௨ భ ்ಽ

ܹ݀ଵ

(28)

8

݀‫ݑ‬ଶ ൌ െ

௨మ

݀‫ݑ‬ଷ ൌ െ

௨య

்ಽ

்ಽ

݀‫ ݐ‬൅ ට

݀‫ ݐ‬൅ ට

തതതതത మ ଶ௨ మ ்ಽ തതതതത మ ଶ௨ య ்ಽ

ܹ݀ଶ

(29)

ܹ݀ଷ

(30)

In these equations ui denotes the fluctuating velocity and TL is the time scale. This model is used in the simulation of particle deposition in turbulent duct flows which resulted in satisfactory predictions by Tian and Ahmadi (2007). Dehbi (2008) used the following Langevin equation to account for the inhomogeneous flow effects: ௨

௨

ఙ

ఙ ்ಽ

݀ ቀ ቁ ൌ െሺ ሻ

ଵ

ଶ

డఙ

்ಽ

డ௫

݀‫ ݐ‬൅ ට ܹ݀ଶ ൅

݀‫ݐ‬

(31)

The last term in Eq. (31) is called drift correction and it is added to avoid the spurious accumulation of fluid particles in region of low turbulent kinetic energy. Dehbi (2008) applied similar formulation based on an idea initially developed by Wilson et al. (1981) to evaluate dispersion of spherical particles for various wall bounded flows. Following Dehbi (2008), the flow field can be divided into two distinct regions, that is the isotropic bulk region and near wall inhomogeneous and anisotropic region. Different formulations for fluctuating velocities are applied in these regions. The following normalized Langevin equation is used near the wall (Dehbi, 2008): ௨

௨

ఙ೔

ఙ೔ ்ಽ೔

݀ ቀ ೔ ቁ ൌ െሺ ೔ ሻ

ଵ

݀‫ ݐ‬൅ ට

ଶ ்ಽ೔

ܹ݀௜ ൅

డఙమ

ߜ ݀‫ݐ‬ డ௫మ ௜ଶ

(32)

whereߜ௜ଶ is the Kronecker delta, andܶ௅௜ denotes the time scales in different directions that will be described subsequently. The normalized Langevin equation given by (32) allows for spatial variation of ߪ௜ and ܶ௅௜ . Iliopoulos and Hanratty (1999) noted that the time scale ܶ௅௜ is typically defined for homogeneous turbulent flows. For inhomogeneous turbulent flows, the accurate estimation of ܶ௅௜ is difficult. In the present study, the expression suggested by Kallio and Reeks (1989) is implemented. That is, the time scale in wall units is given as, ܶ௅ ା ൌ ͳͲǡ݊ା ൑ ͷ

(33) ଶ

ܶ௅ ା ൌ ͹Ǥͳʹʹ ൅ ͲǤͷ͹͵ͳ݊ା െ ͲǤͲͲͳʹͻ݊ା ǡ ͷ ൑  ݊ା ൑ ͳͲͲ and



9

(34)

ఔ

ܶ௅ ൌ ܶ௅ ା ሺ௨‫ כ‬ሻమ

(35)

In the above equations, ݊ା is the wall distance in wall units. That is, ݊ା ൌ

௡௨‫כ‬

(35)

ఔ

ఛ

where‫ כݑ‬ൌ ට is the friction velocity, ߬ is the wall shear stress and ߩ is the fluid density. ఘ For distance close to the wall with ݊ା ൏ ͳͲͲǡ the expressions for the root mean square (RMS) velocity fluctuations as obtained from DNS of turbulent channel flow by Dreeben and Pope (1997) are used. There are, ߪ ାଵ ൌ

ఙభ

ൌ

ߪ ାଶ ൌ

ఙమ

ൌ

ߪ ାଷ ൌ

ఙయ

ൌ

ߝൌ

௨‫כ‬

௨‫כ‬

௨‫כ‬

଴Ǥସ௡శ

(36)

ଵା଴Ǥ଴ଶଷଽሺ௡శ ሻభǤరవల ଴Ǥ଴ଵଵ଺ሺ௡శ ሻమ

(37)

ଵା଴Ǥଶ଴ଷ௡శ ା଴Ǥ଴଴ଵସሺ௡శ ሻమǤరమభ ଴Ǥଵଽ௡శ

(38)

ଵା଴Ǥ଴ଷ଺ሺ௡శ ሻభǤయమమ ଵ

(39)

ସǤହଶଽା଴Ǥ଴ଵଵ଺ሺ௡శ ሻభǤళఱ ା଴Ǥ଻଺଼ሺ௡శ ሻబǤఱ

Here ߪ ାଵ is the nondimensional RMS velocity and ߝ is the dissipation of turbulence kinetic energy. The most critical feature of near wall fluctuation is the quadratic variation of turbulence fluctuation normal to the wall as given by Equation (37). As pointed out by Li and Ahmadi (1993), He and Ahmadi (1999) and Tian and Ahmadi (2007), accounting for quadratic variation of turbulence fluctuation with distance from the wall is critical for proper evaluation of particle deposition rate. Substituting the expression (36)-(38) for the RMS velocities into Eq. (32), the stochastic simulation of fluctuating velocities in the wall layer can be performed. Ghahremani et al. (2014) used the CRW model to investigate the deposition of spherical particles in complex passages of human upper airways. It should be noted that while Eqs. (36)-(38) are developed for the turbulent channel flow, it is applicable to the near wall turbulent flows in other configurations such as pipe flows. As stated before the hydrodynamic torques given by Eqs. (20)-(22) contain velocity gradient components. In order to present a stochastic equation for fluctuating velocity gradient, Kvasnak et al. (2002) used the single point joint velocity-velocity gradient pdf equation. The closed pdf 

10

equations were developed by Ahmadi and Hayday (1988), with the associated Langevin equations for the fluid velocity and velocity gradients. In attempt to present a Lagrangian model for velocity gradient tensor, along the line of Girimaji and Pope (1990), Kvasnak et al. (2002) suggested the following Langevin equation for velocity gradient: ݄݀௜௝ ൌ െ‫ܯ‬௜௝ ݀‫ ݐ‬ᇱ ൅ ‫ܦ‬௜௝௞௟ ܹ݀௞௟ where ݄௜௝ ൌ

డ௨೔ డ௫ೕ

(40)

is the instantaneous velocity gradient tensor.

In Eq. (40) ‫ܯ‬௜௝ is the drift

termǡ ‫ܦ‬௜௝௞௟ is the diffusion coefficient tensor, and ܹ݀௞௟ denotes a tensor-valued Wiener process with zero mean and unit variance.The drift coefficient is expressed as: ‫ܯ‬௜௝ ൌ

௛೔ೕ ିு೔ೕ

(41)

ఛആ

Here ߬ఎ ൌ ሺߥȀߝሻଵȀଶ is the Kolmogorov time scale,‫ܪ‬௜௝ denotes the mean velocity gradient and ߝ is dissipation of turbulence kinetic energy. According to Girimaji and Pope (1990) the diffusion coefficient can be expressed as: ଵ

‫ܦ‬௜௝௞௟ ൌ ܽොඥ݄௠௡ ݄௠௡ ቀߜ௜௞ ߜ௝௟ െ ߜ௜௝ ߜ௞௟ ቁ

(42)

ଷ

It should be noted that Girimaji and Pope (1990) suggested another equation for the components of velocity gradient tensor with different formulation for the drift term. That is, ݀݃௜௝ ൌ ൫െܰ௜௝ ൅ ‫ܯ‬௜௝ ൯݀‫ ݐ‬ᇱ ൅ ‫ܦ‬௜௝௞௟ ܹ݀௞௟

(43)

Here ݃௜௝ is the dimensionless velocity gradient given as, ݃௜௝ ൌ ݄௜௝ ߬ఎ

(44)

where ߬ఎ is the Kolmogorov time scale. In Eq. (43) the drift term contains a deterministic part ܰ௜௝ which can be expressed as: ଵ

ܰ௜௝ ൌ ݃௜௞ ݃௞௝ െ ݃௠௡ ݃௡௠ ߜ௜௝

(45)

ଷ

Here, ‫ ݐ‬ᇱ ൌ

௧ ఛആ

is the dimensionless time and the diffusion coefficient is defined in Eq. (42) where

h( )( ) is replaced with g( )( ). Using stochastic equations (32) and (40) for evaluating the fluctuating velocity and velocity gradient components and substituting these into the translational and rotational particle equations of motion, Eqs. (5)-(8), it is possible to evaluate the trajectory of fibrous particles in the turbulent flows. 

11

2.3. Synthesizing isotropic-homogeneous turbulence To evaluate the accuracy of stochastic modeling for dispersion of ellipsoidal fibers, first the performance of stochastic models is evaluated in the case of an isotropic homogeneous turbulent field. The method proposed by Kraichnan (1970) has been often used to generate Gaussian homogeneous isotropic pseudo turbulent fluctuating fields. Another method for generating fluctuating velocity field is presented by Billson et al. (2004). This method offers a convenient way to prescribe turbulent length and time scales independently. In this approach, a turbulent velocity field is simulated using random Fourier modes. That is, ‫ݑ‬௜ ᇱ ൫‫ݔ‬௝ ൯ ൌ ʹ σே ො ௡ …‘•ሺߢ௝ ௡ ‫ݔ‬௝ ൅ ߰ ௡ ሻߪ௜ ௡ ௡ୀଵ ‫ݑ‬

(46)

where ‫ݑ‬ො௡ , ߰ ௡ and ߪ௜ ௡ denote, respectively, the amplitude, phase and direction of Fourier mode n. Earlier study of Ounis and Ahmadi (1989) reveals that ܰ ൌ ͳͲͲ leads to appropriate fluctuation velocity fields. The synthesized turbulence is generated and used in the analysis. First for each mode n, random anglesߞ ௡ ,ߣ௡ , ߛ ௡ and random phase ߙ ௡ are created. The probability distribution of random phase and random angles are given as: ଵ

ǡͲ ൏ ߞ ௡ ൏ ʹߨ ‫݌‬ሺߞ ௡ ሻ ൌ ቊଶగ Ͳǡܱ‫݁ݏ݅ݓݎ݄݁ݐ‬

(47)

ଵ

ǡͲ ൏ ߣ௡ ൏ ʹߨ ‫݌‬ሺߣ௡ ሻ ൌ ቊଶగ Ͳǡܱ‫݁ݏ݅ݓݎ݄݁ݐ‬ ௡

ଵ

‫݌‬ሺߛ ሻ ൌ ቊଶ

௡

‫݌‬ሺߙ ሻ ൌ

(48)

‫ߛ݊݅ݏ‬ǡͲ ൏ ߛ ௡ ൏ ߨ

(49)

Ͳǡܱ‫݁ݏ݅ݓݎ݄݁ݐ‬

ଵ

ቊଶగ

ǡͲ ൏ ߙ ௡ ൏ ʹߨ

(50)

Ͳǡܱ‫݁ݏ݅ݓݎ݄݁ݐ‬

Then, the highest and smallest wave numbers are specified. Dividing the wave number space into N modes of equally sizeȟߢ , the randomized components of wave number vector are then computed. It should be noted that continuity requires that the unit vector ߪ௜ ௡ and ߢ௝ ௡ would be



12

orthogonal. In the next step a modified von Karman spectrum is selected. The amplitude of each mode ‫ݑ‬ො௡ is then evaluated from, ‫ݑ‬ො௡ ൌ ሺ‫ܧ‬ሺหߢ௝ ௡ หሻȟߢሻଵȀଶ

(51)

The above equations allow computing the turbulent fluctuation components. The fluctuating velocity field is generated for each time step, which leads to zero time correlation. To create correlation in time, a new fluctuating velocity fields are calculated based on an asymmetric time filter given as, ሺܷ௜ ᇱ ሻ௠ ൌ ܽሺܷ௜ ᇱ ሻ௠ିଵ ൅ ܾሺ‫ݑ‬௜ ᇱ ሻ௠

(52)

where m denotes the number of time step andܽ ൌ ‡š’ሺെȟ‫ݐ‬Ȁܶ௅ ሻ, ܾ ൌ ሺͳ െ ܽଶ ሻ଴Ǥହ . The time correlation of new fluctuating field is ‡š’ሺെȟ‫ݐ‬Ȁܶ௅ ሻ . Using these equations for a and b ensures ଶ

ଶ the equality of the variance of velocity fluctuation, i.e. ‫ܷۃ‬௜ ᇱ ‫ ۄ‬ൌ ‫ݑۃ‬௜ ᇱ ‫ۄ‬.

For velocity gradient tensor, each component is calculated by differentiation of the synthesized velocity field. One of the benefits of this formulation is the ability of evaluating the instantaneous flow velocity vector and velocity gradient tensor at the locations of particles, which is important for Lagrangian particle tracking. This method is used in the present study for generating synthesized turbulence fluctuation field. While Fan and Ahmadi (1995a) used the Kraichnan model for dispersion analysis of fibrous particles, to our knowledge this modified procedure was not used earlier for dispersion analysis of spherical particles or elongated fibers in turbulent flows. 3. Results In order to simulate the dispersion of elongated particles, a numerical code was developed for analyzing fiber trajectories in laminar pipe flow, in a synthesized isotropic homogenous turbulence, and turbulent pipe and duct flows. The code includes two options for generation of fluctuating velocity and velocity gradient components. It can synthesize the pseudo-turbulent field according to the formulation of Billson et al. (2004) or use the stochastic modeling for velocity and velocity gradient components with drift correction. In the following sections, the performance of stochastic modeling for velocity and velocity gradient for dispersion of fibers are assessed. An isotropic homogeneous and inhomogeneous turbulent field is synthesized along the



13

fiber trajectories and dispersion of fibers in such flow fields is evaluated.

The stochastic

modeling with the same statistics is also performed and the results for dispersion of fibers are compared with those from the synthesized field. Finally, dispersion and deposition of fibers in a fully developed turbulent pipe and duct flows are evaluated and results are compared with the experimental measurements.

3.1. Validation of fiber motion in laminar pipe flow To check validity of the developed numerical model a test case for laminar pipe flow similar to the studies of Chen and Yu (1991) and Hogberg et al. (2008) is simulated. Figure 3 shows the pipe geometry and the location of inertial coordinate system. In this case a fully developed flow with the average velocity of 0.12 m/s in a pipe with a radius of 1 mm is considered. Here, the fluid viscosity of ߤ= 1.8×10-5 kg/ms, and density of ߩ=1.225 kg/m3 are assumed. An ellipsoidal fiber with the semi minor axis of a=5m, aspect ratio of ߚ=20 and density of ߩ௉ = 1000 kg/m3 is injected into the pipe. The pipe axis is assumed to be aligned with the Y axis of the inertial coordinate system and the gravity is in Z direction. The injection location is at X0= -0.75R and Y0=Z0=0, where R is the pipe radius. The initial orientation of fiber is ߮଴ =45°,ߠ଴ =-45° and ߰଴ =0. The results of numerical simulation for the motion of a single fiber in laminar pipe flow are shown in Fig.4. This figure shows good agreement between the present numerical simulations shown in Fig.4b and the earlier reported results of Hogberg et al. (2008) reproduced in Fig.4a. It should be noted that in the study of Chen and Yu (1991) and Hogberg et al. (2008) the inertial forces and torques were neglected; therefore, some small differences in the trajectory of fiber especially in the Y-Z plane is expected.

The accuracy of the developed computational model for predicting the deposition efficiency of fibers in laminar pipe flows is assessed in this section by comparing the result with the recent study of Feng and Kleinstreuer (2013). In this earlier work, the deposition efficiency of fibers with the equivalent volume diameter of dev=2.41m and various aspect ratios in a pipe with a diameter of d=0.004m and length of L=0.7m for different Reynolds numbers in laminar flow regime were evaluated. Here, the deposition efficiency is defined as the ratio of the number of deposited fibers to the total number of released fibers. The reference coordinate system is shown in Fig.3, gravity is in the Z direction and the lift force is not included in the fiber equations of 

14

motion. For identical conditions, the present computational model is used and the trajectories of 10,000 fibers that are released at the pipe inlet plane are analyzed and the corresponding deposition efficiencies are evaluated for different flow Reynolds numbers. Figure 5 compares the predicted fiber deposition efficiencies with the earlier results of Feng and Kleinstreuer (2013). It is seen that there is a good agreement between the deposition efficiency predicted by the present model and the study of Feng and Kleinstreuer (2013).The slight discrepancy between the present results and the study of Feng and Kleinstreuer (2013) is perhaps due to slight difference in the velocity profile and equation discretization scheme used in the simulations. 3.2. Time evolution of variance of fiber displacement in an isotropic homogeneous field In this section, comparisons of fiber mean-square displacement in the synthesized random velocity field are presented. In these simulations, fibers are injected from a point source atܺ଴ ൌ Ͳ, ܻ଴ ൌ Ͳ, ܼ଴ ൌ Ͳand they are tracked in the isotropic homogeneous synthetic turbulent flow fields. The fibers are released with random orientation and the gravity force is neglected. The random orientation for each fiber is generated using random Euler angles. Each random Euler angle is extracted from a series of random numbers with uniform distribution. Simulations results for the mean and variance of displacements are evaluated and the results are shown in Fig. 6. The results obtained from stochastic modeling using equations (32) and (40) are also evaluated and they are shown in this figure for comparison. It is seen that the stochastic model predictions are in good agreement with the simulated results for the mean and mean-square fiber displacement in the synthetic turbulent flow field. Small discrepancy observed in Fig. 6b is attributed to the statistical error for finite number of fibers. As shown in Fig. 6a, the mean value of displacements are zero because the mean value of velocity components are zero.

3.3. Time evolution of moments of fiber displacement in a fully developed turbulent flow field In this section the mean and variance of translational displacement of fibers in an inhomogeneous turbulent flow field near the wall are evaluated and compared. Fibers are released from a point source at a region near the wall in a fully developed pipe flow and for n+<100 at଴ ൌ ͲǤͲʹͶǡ ଴ ൌ Ͳ, ଴ ൌ Ͳ. For the mean velocity field the logarithmic profile with 

15

—‫ כ‬ൌ ͲǤͶȀ• is used. For simulating the turbulence fluctuations, the synthesized flow field and stochastic modeling are used. The synthesized field produces isotropic instantaneous fluctuation field. Here to account for the anisotropy near wall turbulent flows, the root-mean-square fluctuations of Dreeben and Pope (1997) that are given by equations (36)-(38) are used in the analysis. To obtain various statistics, an ensemble of 1000 fibers are released from the point source and trajectory of each fiber is followed for 1000 time steps with dt=10-6s. The fibers are also released with random orientation and the gravity is neglected. The simulation results for the synthetic turbulent flow are presented in Fig. 7 and compared with those obtained from stochastic model. This figure shows good agreement of the first and second moments of fiber displacements as predicted by the stochastic modeling using equations (32) and (40) and the simulation results of the synthetic turbulent flow field. As noted before, the fibers are released from a point located at଴ ൌ ͲǤͲʹͶ. According to the Fig. 7, the largest mean value and variance is observed for  displacement. Agreement between mean and variance of translational displacements suggests similarity of the fluctuation field generated from the stochastic model with the synthesized pseudo turbulent flow field. Among the two procedures for generation of fluctuating velocity and velocity gradient components, the required computational time of simulations using the stochastic modeling is about 50% of the simulation using the synthesized field. In the following numerical simulations, the stochastic modeling approach is used for the subsequent simulations.

3.4. Deposition of fibers in fully developed turbulent pipe and duct flows Performance of the model for analyzing dispersion and deposition of fibers in fully developed turbulent pipe flow is analyzed in this section. Earlier, Chen and Ahmadi (1997) and Matida et al. (2000) among others, simulated deposition of spherical particles on the pipe wall in turbulent dispersed flows, where a one-way coupled Lagrangian model was used. In their study, they modified the eddy-interaction (DRW) with near wall correction model for evaluating the velocity fluctuations. Their calculations were done in the cylindrical coordinate system. In the present study, similar simulations in the Cartesian coordinate system are performed for elongated particles, except that the CRW model with drift correction is used for the fluctuating velocities.



16

The axial velocity, given by the law of the wall, and RMS values of fluctuating velocities, fitted to DNS data of Dreeben and Pope (1997), are used to represent the mean and fluctuation velocities. The turbulence dissipation rate, , are evaluated using the DNS data of Dreeben and Pope (1997) and the time scales are computed from Eq. (35) as suggested by Kallio and Reeks (1989). The friction velocity which is used for normalization is calculated similar to the study of Fan and Ahmadi (2000). In order to have uniform distribution of fibers at the inlet section, the inlet plane is divided into number of radial and circumferential segments and the fibers are released from these locations. Different particle loadings of 1000, 4000 and 10000 fibers are released at the inlet plane to test the statistical independence of results for dispersion of fibers. Increasing number of fibers to 20000, it was found that the use of 10000 fibers is sufficient to reach statistical accuracy of about 1%. To evaluate the accuracy of simulation results, the geometry and flow Reynolds number are selected to be identical to the experimental study of Shapiro and Goldenberg (1993) for deposition of glass fibers of different size and aspect ratio in a 2 inch diameter horizontal pipe. In these experiment, glass fibers with density of 2230 kg/m3, diameter of d=0.6-2.5m and length of Lp<50m were dispersed in the airflow in the pipe. Deposition velocities for the floor, the ceiling and the side walls were evaluated by means of cellotapes placed on the different sectors of the test section. Cellotapes of 19 mm width was placed at 0° (at the zenith), 90° at the right hand side and 180° in direction of gravity. Each cellotape covered a 43° equivalent arc in the test section. In the present study, the side wall deposition data of Shapiro and Goldenberg (1993) is simulated assuming zero gravity, g=0. In a selected test case, the flow Reynolds number in the experiment was Re=3×104, which corresponds to a friction velocity of‫ כݑ‬ൌ ͲǤ͸݉Ȁ‫ݏ‬. For this friction velocity, the correlations given by Eqs. (36)-(39) are used to evaluate the RMS velocities and dissipation. The time step for particle trajectory analysis needs to be carefully selected. Li and Ahmadi (1993), He and Ahmadi (1999) and Tian and Ahmadi (2007) suggested the use of 0.1p for nonBrownian particles, where,߬௣ ൌ

ఘௗ೛ మ ଵ଼ఓ

is the particle relaxation time and ݀௣ is the diameter of

volume equivalent spherical particle. Dehbi (2008) used min (0.1p, 10-6) for the time step for the numerical simulation of spherical particle motions. For fluctuating velocity gradient 

17

components the scaling time is the Kolmogorov time scale . Girimaji and Pope (1990) used dt=0.1 in the stochastic modeling of velocity gradient components in an isotropic homogeneous turbulence. In the current study, at first a time step of dt=10-5s is selected for numerical simulation. Using this time step, the stochastic modeling of turbulence fluctuations and particle trajectory analysis are performed. For dispersion of larger fibers with relaxation time of 10-5s, the results reveal that decreasing the time step to dt=10-6s, did not affect the dispersion statistics to a noticeable extent. However, for evaluating the deposition velocity of fibers with lower aspect ratios, using a time step of dt=10-6s leads to a more accurate model prediction in comparison with the experimental data. Therefore in order to satisfy all of the mentioned constraints on the time step size the deposition data are obtained using the time step of min (0.1p, 10-6). The deposition mechanism of fibers is more complicated than spherical particles due to its elongated shape. In the current study, the condition for deposition proposed by Fan and Ahmadi (1995b) is used. Accordingly, if the distance between the fiber mass center and pipe wall is less than the fiber semi minor axis length, a, the fiber will deposit on the wall. If the distance between the mass center of particle and the pipe wall is larger than the fiber semi major axis length, b, the fiber will not reach the wall. When the distance between the mass center of fiber and the wall is smaller than b and larger than a the deposition of fiber depends on the orientation of major axis.

According to derivation of Fan and Ahmadi (1995b) the possible deposition of an ellipsoidal fiber can be identified by: ೎೚ೞഀ మ

ೞ೔೙ഀ మ

ሺ௕మ ି௔మ ሻమ ௦௜௡మ ఈǤ௖௢௦ మ ఈ൤ቀ ೌ ቁ ାቀ ್ ቁ ൨

ߦൌඨ

భ

భ

ߞൌ

(53)

భ మ

ଵା௔మ ௕మ ቀ మ ି మ ቁ ௦௜௡మ ఈǤ௖௢௦మ ఈ ೌ ್ ೎೚ೞഀ మ

భ

ೞ೔೙ഀ మ

഍ మ

ିకቀ మ ି మ ቁ௦௜௡ఈǤ௖௢௦ఈିටቀ ೌ ቁ ାቀ ್ ቁ ିቀೌ್ቁ ೌ ್ ೎೚ೞഀ మ

(54)

ೞ೔೙ഀ మ

ቀ ೌ ቁ ାቀ ್ ቁ

ߙ ൌ ܿ‫ି ݏ݋‬ଵ ሺ‫ܣ‬ଷଵ ଶ ൅ ‫ܣ‬ଷଷ ଶ ሻ଴Ǥହ

(55)

The coordinate systems and the deposition situation for an ellipsoidal fiber are shown in the Fig.8. A positive value of ߞindicates no contact between the fiber and the wall, while a negative



18

value indicates that the fiber touches the wall and stick the wall. In the present study the effect of particle rebound is neglected. Results of the numerical simulations for dispersion and deposition of fibers are presented in Figs. 9-12. In Fig.9, 1000 fibers are released into the pipe with random orientations at zero translational and rotational velocities and the equations of motion given by (5)-(8), are solved for each fiber until the fiber is either deposited on the pipe wall or leaves the pipe. The fiber positions along their trajectories are stored and plotted at every 100 time steps in Fig.9. This figure shows the simulated dispersion pattern of fibers in a fully developed turbulent pipe flow for the conditions similar to the experimental study of Shapiro and Goldenberg (1993). In particular, traces of fibers moving with the flow are shown in this figure. It is seen that the fibers roughly follow the mean flow but they are dispersed by the turbulence fluctuations and are deposited on the wall. Various approaches for evaluating the deposition velocity in computer simulations were suggested in the literature. The traditional approach as suggested by Li and Ahmadi (1993) and Fan and Ahmadi (2000) is given as, ‫ݒ‬ௗ௘௣ ൌ

௃ ஼బ

ൌ

οே೏ ൗο௧ ೏ ேబ ൗ௛ బ

(56)

where ‫ ܬ‬denotes the particle mass flux to the wall per unit time,‫ܥ‬଴ is the concentration of particles near the surface,οܰௗ is the total number of deposited particles during time interval ο‫ݐ‬ௗ , ܰ଴ is the initial number of released particles which are uniformly distributed within distance ݄଴ near the wall. Here ݄଴ corresponds to the pipe radius in the present study. Matida et al. (2000) applied the following formulation for evaluating the deposition velocity of spherical particles in the turbulent pipe flow: ‫ݒ‬ௗ௘௣ ൌ

ௗ೟ തതതത ௨೑ ସ௅

ி

Žሺ భሻ

(57)

ிమ

In this equation ‫ܨ‬ଵ and ‫ܨ‬ଶ are, respectively, the concentration of particles at the inlet and outlet planes,݀௧ is the pipe diameter,‫ ܮ‬is the pipe length, and ‫ݑ‬ തതത௙ denotes the average fluid velocity in the pipe. This formulation assumes a constant velocity and concentration profiles across the pipe.



19

Here both approaches given by Eqs. (56) and (57) are used to compute the deposition velocity of fibers and the results are compared with the experimental data of Shapiro and Goldenberg (1993) in Fig.10. In this figure Lp is the major axis length of ellipsoidal fiber. The significant increase in the deposition velocity with aspect ratio is observed from this figure. This is due to the increase in the wall capture efficiency due to interception. As mentioned by Mehel et al. (2010), the formulation for deposition velocity given by Eq. (57) overpredicts the deposition velocity for spherical particles. Similar overestimation is observed from Fig. 10 for the fiber deposition especially when smaller number of fibers is used. It is conjectured that the overestimation of deposition velocity is due to statistical error for number of released fibers, which decreases with the increase in the number of fibers. Clearly the accuracy of the computational model prediction for the deposition velocity increases by increasing the number of fibers. As noted before, increasing the number of fibers beyond 10000 did not affect the deposition data considerably. Fig. 10 also shows that the predicted deposition velocities by both Eqs. (56) and (57) for NP=10000 are lower than the experimental data for 40 and 50 μm long fibers. This could be in part due to the experimental uncertainty in this range, or slight difference between the velocity profiles in the present study and experiment of Shapiro and Goldenberg (1993). It should also be pointed out that the data of Shapiro and Goldenberg (1993) are for cylindrical fibers, while the present analyses were for ellipsoidal fibers. In summary, Eq. (56) seems to lead to more accurate estimate of deposition velocity of fibers for aspect ratios smaller than 15.

The floor deposition velocity of glass fibers in a turbulent pipe flow for various flow Reynolds numbers are presented in the study of Shapiro and Goldenberg (1993). The measured deposition velocities are presented for fibers with different equivalent relaxation times, ߬ ା ௘௤ ൌ

௠௨‫כ‬

మ

ଵ

ఓఔ ௄ᇲ

‫ ܭ‬ᇱ ൌ ͵ሺ‫ ܭ‬ᇱ ௫ ᇲ ௫ ᇲ

where‫ ܭ‬ᇱ denotes the orientation averaged resistance:

ିଵ

൅ ‫ ܭ‬ᇱ ௬ᇲ௬ᇲ

ିଵ

൅ ‫ ܭ‬ᇱ ௭ᇲ௭ᇲ

ିଵ ିଵ

ሻ

(58)

where‫ ܭ‬ᇱ ሺሻᇲ ሺሻᇲ are defined in Eqs.(16)-(17). The non-dimensional deposition velocity is defined as the ܸ ା ൌ ‫ݒ‬ௗ௘௣ Ȁ‫ כݑ‬. Numerical simulations are carried out for different cases similar to the condition of experimental study and the corresponding deposition velocities are plotted in Fig.



20

11. For each case, the friction velocity is evaluated based on the flow Reynolds number. For a range of fiber aspect ratios, a series of simulations was performed and the corresponding deposition velocity is computed using Eq. (56). Fig. 11 shows good agreement between the numerical simulation results and the experimental data of Shapiro and Goldenberg (1993). For Re=30000, the deposition velocity increases sharply with ߬ ା ௘௤ . As the flow Reynolds number increases, the rate of increase of deposition velocity with the equivalent relaxation time decreases. For Re=87000, the deposition velocity of glass fibers, in the studied size range, is approximately constant,ܸ ା ൌ ͲǤͳͳ. That is, for this size range the deposition velocity reaches to the constant value independent of the increase in fiber length. The constant non-dimensional deposition velocity for the higher Reynolds number and thicker fibers can be explained by noting that the thicker fibers have sufficiently high inertia to penetrate through the viscous sublayer and deposit on the wall and further increase in the particle inertia does not affect the deposition process. It should be also noted that a constant ܸ ା means the deposition velocity increases linearly with increase of shear velocity, u*. Comparison between the order of magnitude for different forces at the Reynolds number of Re=87000 revealed that the drag force (mean flow drag) is about 100 times greater than the gravity force; therefore, most fibers are carries out to the pipe exit by the flow rather than depositing on the wall. The constant trend of nondimensional deposition velocity for spherical particles with high inertial was observed in the experimental study of Liu and Agarval (1974) and the numerical simulation of Fan and Ahmadi (1995a) and Tian and Ahmadi (2007). Fan and Ahmadi (2000) and Zhang et al. (2001) observed similar constant trends for high inertia ellipsoidal fibers. In Fig.11 simulation results for non-dimensional deposition velocity in the absence of turbulent velocity gradient fluctuation are also presented for flow Reynolds numbers of 30000, 54000 and 87000 for comparison. It is seen that the effect of fluctuating velocity gradient is quite small for the range of simulated flow condition and fibers sizes. It may then be concluded that slightly more accurate results for deposition velocity are obtained when the effect of random velocity gradient is included in the computational model, but neglecting the velocity gradient fluctuation will not have a significant effects. In this section, the accuracy of the developed computational model for predicting the dispersion and deposition of fibers in turbulent duct flows is assessed by comparison of the results with the experimental data of Kvasnak and Ahmadi (1995) for deposition rate of glass



21

and paper fibers. The wind tunnel used in the experiment had a test section with cross sectional area of 15.25×2.54 cm2 that was 36 cm long. Kvasnak and Ahmadi (1995) also suggested the following empirical equation for non-dimensional deposition velocity (‫ݑ‬ା ௗ ൌ

௨೏ ௨‫כ‬

) on the duct

floor for both spherical and elongated particles, ା ଶ ା ‫ݑ‬ା ௗ ൌ ͶǤͷ ൈ ͳͲିସ ሺ߬௩௢௟ ሻ ൅ ߬௩௢௟ ݃ା ൅ ͷ ൈ ͳͲିଷ ሺ‫ܮ‬ା ሻଶ

Here ‫ܮ‬ା ൌ ‫ כݑܮ‬Ȁߥ is the dimensionless particle length, ݃ା ൌ

(59) ௚ఔ

௨‫ כ‬య

is the non-dimensional gravity,

ା ݃ is the gravitational acceleration, ߥ is the fluid kinematic viscosity, and ߬௩௢௟ denotes the

equivalent volume sphere relaxation time that for an ellipsoid of revolution is expressed as, ା ߬௩௢௟ ൌ

మ

ௌௗ శ ఉమȀయ

(60)

ଵ଼

Here S is particle to fluid density ratio, d=2a and݀ ା ൌ ݀‫ כݑ‬Ȁߥ. In Eq. (59) the first term on the right-hand side is the turbulent eddy-impaction effects, the second term is the gravitational sedimentation and the third term is attributed to the enhanced deposition due to fiber elongated shape. In the present study, for the case identical to the experimental study of Kvasnak and Ahmadi (1995), the deposition velocity is simulated and results are presented in Fig.12 with the experimental data.

For fibers with a given diameter, both computer simulations and

experimental data show sharp increase of deposition velocity with an increase of fiber length. That is, as the fiber aspect ratio increases, the wall capture efficiency due to interception increases. Fig. 12 shows reasonable agreement of the present simulation results with the experimental data of Kvasnak and Ahmadi (1995) and the empirical equation prediction for fibers with L<50m. For longer fibers with L>50m, however, the simulation results underpredict the experimental deposition velocity. This may be attributed to difference between velocity profiles in the experiment and numerical simulation. Also it should be noted that in the present study the creeping flow formulations for hydrodynamic forces and torques are used. However, for long fibers with L>50m, the flow condition may be non-creeping and this may affect the deposition results. In Fig. 12, the simulation results for deposition velocity in the absence of fluctuation velocity gradient are also shown for comparison. It is seen that the simulations without accounting for the stochastic variation of velocity gradient slightly overpredict the deposition velocity for fibers



22

with lower aspect ratio. Neglecting the velocity gradient fluctuation for fibers with higher aspect ratio, however, leads to underprediction of deposition velocity. The presented comparisons of the computational model predictions with the experimental data show that the developed numerical method can properly predict the deposition of fibers in the turbulent pipe and duct flows. The importance of fluctuating velocity gradient on the dispersion and deposition of fibers, however, requires further investigation for other geometries and other flow conditions.

3.5 Effect of lift force A series of simulations were performed to assess the importance of lift force. In the study of Fan and Ahmadi (1995) it was reported that the lift force is important for deposition of fibers in the vertical duct flows. Accordingly, for upward and downward turbulent pipe and duct flows, the deposition of fibers is controlled by the velocity fluctuations and is noticeably affected by the lift force. In the present study, however, the numerical simulations are performed for the horizontal pipe flows and the influence of lift force on deposition is rather difficult to interpret. For floor deposition of fibers in horizontal pipe flows, the lift force may increase or decrease the deposition depending on the direction of relative velocity,ሺ࢛ െ ࢜ሻ. The results of the present study on the influence of the lift force on deposition of fibers in the turbulent pipe flows are shown in Figs. 10-12. Fig.10 shows that the influence of the lift force on the deposition velocity of fibers on the side walls is not noticeable. Inspecting Eq. (5), it is clear that the magnitude of the lift force depends on the size of fibers as well as the shear rate. For the condition of Fig.10, the shear rate and size of the fibers are not sufficiently high to produce a significant lift force; therefore, roughly the same deposition velocities are obtained when the lift force is included or neglected in the analysis. For floor deposition of fibers, numerical simulations are performed for higher Reynolds numbers and for fibers with higher aspect ratios. Fig. 11 shows that the influence of the lift force on the deposition velocity of fibers increases for higher Reynolds numbers corresponding to the higher shear rates. Nevertheless, it is seen that the computed deposition velocities without the lift force are comparable with those with the lift force included in the analysis. For the case of turbulent duct flow shown in Fig.12, the floor deposition velocity for longer fibers increases slightly with the inclusion of the lift force.



23

4. Conclusions In this study a one-way coupling Eulerian-Lagrangian computational method for analyzing dispersion and deposition of fibers in laminar and turbulent duct flows was developed. For the case of turbulent flow, the influences of turbulent fluctuations were included using appropriate stochastic models for the velocity gradients in addition to the fluctuating velocities. To generate turbulence fluctuations, two different methods were used. In the first method, synthetic continuous isotropic homogenous, as well as, inhomogeneous Gaussian random fields were simulated and the corresponding instantaneous velocities and velocity gradients were evaluated. In the second procedure, a stochastic filtered white noise model was used for generating instantaneous velocity and velocity gradient components. The second approach offers lower computational time. For inhomogeneous near wall region, the root-mean square fluctuating velocities and timescales were obtained from the fitting to the available DNS data of turbulent channel flow. A series of simulations were performed and fluctuating velocity and velocity gradient fields were generated and used for fiber transport and deposition analysis. The computation model was validated by comparison of the model predictions with the available experimental data. Based on the presented results the following conclusions are made: 1. The simulation procedure developed for dispersion analysis of fibers in laminar and turbulent flows seems to lead to reasonable results. 2. Using the synthesized field and stochastic modeling for generation of fluctuating velocity components, similar dispersion patterns are obtained in isotropic homogenous turbulence and inhomogeneous turbulent pipe flows. 3. The calculated deposition velocities for fibers with ߚ ൑ ͳͷin fully developed pipe and channel flows are in good agreement with the experimental data. The increase of fiber aspect ratio increases the deposition velocity for both floor and side deposition in turbulent pipe and channel flows. 4. At Re=87000, constant deposition velocity is obtained for floor deposition in turbulent pipe flow for ellipsoidal fibers with ߬ ା ௘௤ ൐ ͺ.



24

5. The influence of fluctuating velocity gradient on the fiber deposition velocity is found to be small for the range of parameters studied. However, slightly more accurate deposition velocities are obtained by including the fluctuation velocity gradient in the analysis. 6. The influence of the lift force on the deposition velocity of fibers on the side walls is negligibly small. For the floor depositions, however, the lift force slightly increases the deposition velocity for longer fibers.

Acknowledgement The third author appreciates support from Iran’s National Elite Foundation.

References Ahmadi, G. (1970). Analytical prediction of turbulent dispersion of finite size particles. Ph.D. Thesis, Purdue University. Ahmadi, G., & Hayday, A. A. (1988). A probability density closure model for turbulence. Acta Mechanica, 72, 55-71. Billson, M., Eriksson, L. E., & Davidson, L. (2004). Modeling of synthetic anisotropic turbulence and its sound emission. 10thAIAA/CEAS Aeroacoustics Conference, AIAA, Manchester, United Kingdom. Brenner, H. (1963). The Stokes resistance of an arbitrary particle. Chemical Engineering Science, 18, 1–25. Brooke, J. W., Kontomaris, K., Hanratty, T. J., & McLaughlin, J.B. (1992). Turbulent deposition and trapping of aerosols at a wall. Physics of Fluids A, 4, 825-834. Chen, Q., & Ahmadi, G. (1997). Deposition of particles in a turbulent pipe flow. Journal of Aerosol Science, 28, 789-796. Chen, Y., &Yu, C. (1991). Sedimentation of fibers from laminar flows in a horizontal circular duct. Aerosol Science and Technology, 14, 343-347. Cleaver, J. W., & Yates, B. (1975). A Sublayer model for deposition of particles from turbulent flows. Chemical Engineering Science, 30, 983-992. Dastan, A., & Abouali, O. (2013). Microfiber motion and web formation in a microchannel heat sink: a numerical approach. Computer and Fluids, 71, 28-40.



25

Dastan, A., Abouali, O., & Ahmadi, G. (2014). CFD simulation of total and regional fiber deposition in human nasal cavities. Journal of Aerosol Science, 69, 132-149. Dehbi, A., (2008). Turbulent particle dispersion in arbitrary wall-bounded geometries: a coupled CFD-Langevin-equation based approach. International Journal of Multiphase flow, 34, 819-828. Dreeben, T. D., & Pope, S. B. (1997). Probability density function and Reynolds-stress modelling of near-wall turbulent flows. Physics of Fluids, 9, 154-163. Fan, F. G., & Ahmadi, G. (1993). A sublayer model for turbulent deposition of particles in vertical ducts with smooth and rough surfaces. Journal of Aerosol Science, 24, 45-64. Fan, F. G., & Ahmadi, G. (1994). On the sublayer model for turbulent deposition of particles in presence of gravity and electric fields. Aerosol Science and Technology, 21, 49-71. Fan, F. G., & Ahmadi, G. (1995a). Dispersion of ellipsoidal particles in an isotropic pseudoturbulent flow field. Journal of Fluids Engineering, 117, 154-161. Fan, F. G., & Ahmadi, G. (1995b). A sublayer model for wall deposition of ellipsoidal particles in turbulent stream. Journal of Aerosol Science, 26,813-840. Fan, F., & Ahmadi, G. (2000). Wall deposition of small ellipsoids from turbulent air flow - A Brownian dynamics simulation. Journal of Aerosol Science, 31, 1205-1229. Feng, Y., & Kleinstreuer, C. (2013). Analysis of non-spherical particle transport in complex internal shear flows. Physics of Fluids, 25, 091904. Fichman, M., Gutfinger, C., & Pnueli, D. (1988). A Model for Turbulent Deposition of Aerosols. Journal of Aerosol Science, 19, 123-136. Friedlander, S. K., & Johnstone, H. F. (1957). Deposition of suspended particles from turbulentgas streams. Industrial and Engineering Chemistry, 49, 1151. Gallily, I., & Cohen, A. (1979). On the orderly nature of the motion of nonspherical aerosol particles. II. Inertial collision between a spherical large droplet and an axially symmetrical elongated particle. Journal of Colloid and Interface Science, 68, 338-356. Girimaji, S. S., & Pope, S. B. (1989). A diffusion model for velocity gradient in turbulence. Physics of Fluids A, 2, 242-256. Ghahremani, E., Abouali, O., Emdad, H., & Ahmadi, G. (2014). Numerical analysis of stochastic dispersion of micro-particles in turbulent flows in a realistic model of human nasal/upper airway. Journal of Aerosol Science, 67, 188–206. Goldstein, H. (1980). Classical mechanics, 2nd edition. Reading, Addison-Wasley Company: MA, USA.



26

Gosman, A. D., & Ioannides, E. (1983). Aspect of computer simulation of liquid-fuelled combustors. Journal of Energy, 7, 482-490. Harper, E.Y., & Chang, I. D. (1968). Maximum dissipation resulting from lift in a slow viscous shear flow. Journal of Fluid Mechanics, 33, 703-727. He, C., & Ahmadi, G. (1999). Particle deposition in a nearly developed turbulent duct flow with electrophoresis. Journal of Aerosol Science, 30, 739-758. Hogberg, S., Akerstedt, H. O., Lundstorm, T. S., Freund. J. B. (2008). Numerical model for fiber transport in the respiratory airways. Proceedings of 19th International symposium on Transport Phenomena. Reykjavic, Iceland. Hutchinson, P., Hewitt, G. F., & Dukler, A. E. (1971). Deposition of liquid or solid dispersions from turbulent gas streams: a stochastic model. Chemical Engineering Science, 26. 419-439. Iliopoulos, I., & Hanratty, T. J. (1999). Turbulent dispersion in non-homogeneous field. Journal of Fluid Mechanics, 293, 45-71. Jeffery, G. B. (1922). The Motion of ellipsoidal particles immersed in a viscous fluid. Proceedings of the Royal Society of London. Series A, 102, 161-179. Kallio, G. A., & Reeks, M. W. (1989). A numerical simulation of particle deposition in turbulent boundary layers. International Journal of Multiphase flow, 15, 433-446. Kleinstreuer, C., & Feng, Y. (2013). Computational analysis of non-spherical particle transport and deposition in shear flow. Journal of Biomechanical Engineering, 135, 002108. Kraichnan, R. H., (1970). Diffusion by random velocity field. Physics of Fluids, 11, 22-31. Kvasnak, W., & Ahmadi, G. (1995). Fibrous particle deposition in a turbulent channel flow-An experimental study. Aerosol Science and Technology, 26, 641-652. Kvasnak, W., Ahmadi, G., & Schmidt. D. J. (2002). An engineering model for fuel spray formation of deforming droplets. Atomization and Sprays, 14, 289-339. Li, A., & Ahmadi, G. (1993). Deposition of aerosol on surfaces in a turbulent channel flow. International Journal of Engineering Science, 31, 435-451. Liu, B. Y. H., & Agarval, J. K. (1974). Experimental observation of aerosol deposition in turbulent flow. Journal of Aerosol Science, 5, 145-155. Matida, E. A., Nishino, K., & Torii, K. (2000). Statistical simulation of particle deposition on the wall from turbulent dispersed pipe flow. International Journal of Heat and Fluid Flow, 21, 389402. Maxey, M. R. (1987). The gravitational settling of aerosol particle in homogeneous turbulence and random flow fields. Journal of Fluid Mechanics, 174, 441,445.



27

Mclaughlin, J. B. (1989). Aerosol Particle Deposition in Numerically Simulated Channel Flow. Physics of Fluids A, 7, 1211-1224.

Mehel, A., Taniere, A., Oesterle, B., & Fontaine, J.R. (2010). The influence of an isotropic Langevin dispersion model on the prediction of micro- and nanoparticle deposition in wallbounded turbulent flows. Journal of Aerosol Science, 41, 729-744. Ounis, H., & Ahmadi, G. (1989). Motions of small rigid spheres in a simulated random velocity field. Journal of Engineering Mechanics, 115, 2107-2121. Ounis, H., & Ahmadi, G. (1990). Analysis of dispersion of small spherical particles in a random velocity field. Journal of Fluids Engineering, 112, 114-120. Ounis, H., & Ahmadi, G. (1991). Motions of small particles in simple shear flow field under microgravity condition. Physics of Fluids A, 3, 2559-2570. Reeks, M. W., (1977). On the dispersion of small particles suspended in an isotropic turbulent flow. Journal of Fluid Mechanics, 83, 529-546. Shapiro, M., & Goldenberg, M. (1993). Deposition of glass fiber particles from turbulent air flow in a pipe. Journal of Aerosol Science, 24, 65-87. Soltani, M., & Ahmadi, G. (1995). Direct numerical simulation of particle entrainment in turbulent channel flow. Physics of Fluid A, 7, 647-657. Soltani, M., & Ahmadi, G. (2000). Direct numerical simulation of curly fibers in turbulent channel flow. Aerosol Science and Technology, 33, 392-418. Soltani, M., Ahmadi, G., Ounis, H., & McLaughlin, J. B., (1998). Direct numerical simulation of charged particle deposition in a turbulent channel flow. International Journal of Multiphase flow, 24, 77-92. Squires, K. D., & Eaton, J. K. (1991). On the preferential concentration of solid particles in turbulent. Physics of Fluids A, 3, 1169-1178. Tian, L., & Ahmadi, G. (2007). Particle deposition in turbulent duct flows- comparisons of different model predictions. Journal of Aerosol Science, 38, 377-397. Tian, L., Ahmadi, G., Wang, Z., & Hopke, P. (2012). Transport and deposition of ellipsoidal fibers in low Reynolds number flows. Journal of Aerosol Science, 45, 1-18. Tian, L., & Ahmadi, G. (2013). Fiber transport and deposition in human upper tracheobronchial airways. Journal of Aerosol Science, 60, 1-20. Wilson, J. D., Thurtell, G. W., & Kidd, G. E. (1981). Numerical simulations of particle trajectories in inhomogeneous turbulence. Part2: systems with variable turbulent velocity scale. Boundary Layer Meteorology, 27, 163-169.



28

Zhang, H., & Ahmadi, G. (2000). Aerosol particle transport and deposition in vertical and horizontal turbulent duct flows. Journal of Fluid Mechanics, 406, 55-80. Zhang, H., Ahmadi, G., Fan, F., & McLaughlin, J. B. (2001). Ellipsoidal particles transport and deposition in turbulent channel flows. International Journal of Multiphase flow, 27, 971-1009. Zhang, H., Ahmadi, G., & Asgharian, B. (2007). Transport and deposition of angular fibers in turbulent channel flows. Aerosol Science Technology, 41, 529-548.

Figure captions: Fig. 1. Three coordinate systems for description of a fiber motion Fig. 2. Geometry of an ellipsoidal fiber Fig. 3. Pipe geometry and reference coordinate system Fig. 4. Comparison of computed fiber trajectories with the study of Hogberg et al (2008) a) Reproduced from Hogberg et al. (2008), b) Present study Fig. 5. Comparison between deposition efficiency of fibers in laminar pipe flow Fig. 6. Comparison of fiber mean-square translational displacements in anisotropic homogeneous turbulent flow for, NP=1000, a=2m, b=c=50m,‫ כݑ‬ൌ ͲǤͶ݉Ȁ‫ݏ‬, initial location, ܺ଴ ൌ Ͳ, ܻ଴ ൌ Ͳ, ܼ଴ ൌ Ͳ a) Mean value of displacement, b) Variance of displacement Fig. 7. Comparison of first and second moments of translational displacement in a fully developed pipe flow for, NP=1000, a=2m, b=c=50m released from initial location ofܺ଴ ൌ ͲǤͲʹͶ݉, ܻ଴ ൌ Ͳ, ܼ଴ ൌ Ͳ for ‫ כݑ‬ൌ ͲǤͶ݉Ȁ‫ݏ‬, R=0.025m a) Mean value of displacement, b) Variance of displacement Fig. 8. Fiber deposition scheme on a smooth surface Fig. 9. Dispersion pattern of fibers in turbulent pipe flow for, Re=30000, NP=1000, u*=0.6, a=1m, dt=10-5s, =20, ߬ ା ௘௤ ൌ ͻǤ͹͸, Red. Dispersion pattern of fibers if moving with mean flow, Black



29

Fig. 10. Comparison of the simulated side deposition of glass fibers with different length, LP in a turbulent pipe flow for, Re=30000, u*=0.6m/s, a=1.9m, dt=10-6s, and with initial random orientation Fig. 11. Comparison of floor deposition velocity of glass fibers in a turbulent pipe flow, a=1.5m, initial random orientation Fig. 12. Comparison of floor deposition of glass fibers with different length in a turbulent duct flow for Re=17000, a=2.5m, dt=10-6s, initial random orientation



30

Fig.1:

  



             



31

Fig. 2:

                     



32

Fig. 3:

  



          





33

Fig. 4a:

Fig. 4b:

        





34

Fig. 5: 

   





35

Fig. 6a:

         



36

Fig. 6b:



37

Fig. 7a:





38

Fig.7.b:

         



39

Fig.8:

   



            



40

Fig.9:

     



     



41

Fig.10:

        



  

42

Fig.11:

            

43

Fig.12:

 

Highlights x x x x x



Dispersion and deposition of ellipsoidal fibers are evaluated in turbulent flows. Lagrangian approach is used to simulate the motion of fibers. The influence of turbulent fluctuations is incorporated using stochastic modeling. Performance of stochastic models for turbulent fluctuations is evaluated.

Influence of fluctuating velocity gradient on the fiber deposition velocity is small.

44