Comparison of a new Eulerian model with a modified Lagrangian approach for particle distribution and deposition indoors.

Understanding of aerosol dispersion characteristics has many scientific and engineering applications. It is recognized that Eulerian or Lagrangian approach has its own merits and limitations. A new Eulerian model has been developed and it adopts a simplified drift-flux methodology in which external forces can be incorporated straightforwardly. A new near-wall treatment is applied to take into account the anisotropic turbulence for the modified Lagrangian model. In the present work, we present and compare both Eulerian and Lagrangian models to simulate particle dispersion in a small chamber. Results reveal that the standard k-ε Lagrangian model over-predicts particle deposition compared to the present turbulence-corrected Lagrangian approach. Prediction by the Eulerian model agrees well with the modified Lagrangian model.

Introduction

Gaining understanding of particle transport and deposition indoors has numerous engineering applications. There has been considerable interest over the past decade in exposure assessment of fine particles inhalation and its influence on public health. The anthrax mailing accidents following the terrorist attacks of 11 September 2001 have generated enormous concern in design and application of ventilation strategy for protecting indoor environment against the intentional release of biological agents. Phrases like aerobiological engineering or immune-building technology are coined recently to reflex the importance of studying aerosol behaviors indoors. Better understandings of the aerosol dynamics in terms of mixing time and dispersion rate are vital to decide the positioning of air toxic sensors which have became an important element for monitoring buildings (Gadgil et al., 2003). Improper placement of the sensors can impair the ability of the sensors to provide the first response decision. In addition, understanding aerosol dispersion and transport is very essential in the prevention of nosocomial transmission of airborne pathogens (Cole and Cook, 1998; Li et al., 2005). After the outbreak of severe acute respiratory syndrome (SARS) in South East Asia 2003, there is increasing research interest in studying transport and control of airborne bacteria or viruses in indoor environments (Beggs et al., 2005; Nicas et al., 2005) and in confined environment like flight cabinet (Mangili and Gendreau, 2005). While conducting in situ experiment can provide detailed information on microorganisms transport and survival in enclosed environments, potential danger and cost involved should be carefully considered. Computational fluid dynamics (CFD) provides a very cost-effective way to perform parametric studies to investigate the microorganism dynamic behavior prior to full-scale experiments.

There are two modeling approaches for two-phase flow problems; namely the Eulerian–Eulerian model (hereafter refers to Eulerian model) and the Eulerian–Lagrangian model (hereafter refers to Lagrangian model). The Eulerian method considers the particle phase as another continuum. Governing equations derived from the mass (species) conservation condition are solved to give details of the particle concentration field. However, some studies treated the particle phase as scalar species (no inertia) and the results should be inferred cautiously (Lu et al., 2005; Noakes et al., 2004). Gravitational settling was considered in some Eulerian CFD models (Murakami et al., 1992; Holmberg and Li, 1998; Zhao et al., 2004), nevertheless, the deposition rate in those studies was only estimated empirically or simply ignored. It has been shown that, in certain circumstances, ignoring the deposition flux may result in numerical instability problems. Theoretical model needs to be developed to evaluate deposition rate according to local turbulent flow condition.

The second approach is the Lagrangian method, which treats the dynamics of a single particle by the trajectory method. Under this framework, generally speaking, the flow field is obtained by applying Reynolds averaged Navier–Stokes (RANS) turbulent models. The flow and other quantities obtained are ensemble-averaged components. Near-wall fluctuating velocities are highly anisotropic with the component normal to the wall substantially smaller than those in the other two directions. Proper treatment for the fluctuation components is critical to particle deposition modeling. However, many previous RANS approaches for indoor environments ignored this effect and lead to incorrect deposition rates.

Equation of motion resulting from various forces exerting on an individual particle is solved to acquire the single-particle trajectory. A large number of sample particles should be analyzed before statistical conclusions can be drawn. A number of Lagrangian simulations have been carried out for particle transport and deposition in ventilated single-zone room (Zhao et al., 2004), two-zone chamber (Lu et al., 1996) and multi-zone chamber (Chung, 1999). Unfortunately, the effect of turbulence on particle phase was overlooked in most of those models, even though the airflow fields were simulated with turbulence models. Recently, large eddy simulation (LES) has applied to simulate particle transport and deposition for simple single-zone geometries (Béghein et al., 2005; Bouilly et al., 2005).

Recently, the present authors have proposed a new Eulerian model which takes external drift forces into account (Chen et al., 2006). Gravitational settling has been incorporated as an external drift velocity and the results agree well with the experimental measurements. Lately, a new Lagrangian model has been developed to improve turbulent intensity in the vicinity of the wall (Lai and Chen, 2007). In the present work, we compared particle distribution and deposition rates for a small model chamber by the two approaches.

Modeling approaches and geometry

The geometry of the single-zone model room is shown in Fig. 1 . The room dimension are length (x)×width (y)×height (z)=0.8 m×0.4 m×0.4 m. Its inlet and outlet are of the same size, 0.04 m×0.04 m. Their centers are located at x=0, y=0.2 m, z=0.36 m and x=0.8 m, y=0.2 m, z=0.04 m, respectively. The symmetrical plane at y=0.2 m is referred as the center plane in the following discussion. Two inlet velocities, 0.225 and 0.45 m s−1 (corresponding to air exchange rates of 10 and 20 h−1, respectively), are tested. The room air temperature is set as 27 °C. Due to the high computational cost for continuously particle tracking for Lagrangian models, particles were injected only once instead of continuous injection adopted in the Eulerian approach. This is a very common practice for almost all Lagrangian simulations; nevertheless, due to the different nature of particle injection, it imposes a key constraint comparing to the Eulerian models.

Fig. 1

Schematics of the model room.

The main objective of the present work is to highlight and compare the two approaches on the prediction of particle phase dispersion in a chamber, thus the air phase model is not mentioned here and the approach can be found elsewhere (Chen et al., 2006; Lai and Chen, 2007). In brief, the airflow field was resolved by a RNG k–ε turbulence model for both Eulerian or Lagrangian approaches. The simulations were performed with the aid of the commercial CFD code FLUENT (2005).

Eulerian model

A simplified Eulerian drift–flux model has been developed to take full advantage of the extremely low volume fraction of indoor particles. The term “drift–flux” (or drift velocity) stands for particle flux (or velocity) caused by effects other than convection, i.e. gravitational settling and diffusion for the current work. As the convective velocity of the particle phase is the same as the air phase, the complexity of the two-phase flow system is greatly reduced. The governing equation for particle transport in turbulent flow field is given aswhere u is the air phase velocity vector, C is the particle mass concentration, kg m−3 (or number concentration, m−3) of particle size group i (hereafter in this paper, the subscript i denotes particle size group), is the particle settling velocity, ε p is the particle eddy diffusivity, and D is the Brownian diffusion coefficient and is the mass concentration source term.

For coarse particles, the loss by deposition (i.e. sedimentation) must be properly treated. In the current approach, the concentration field is divided into two regions: the core region and the concentration boundary layer. The methodology adopted here is to get the distribution of particles in the core region with the three-dimensional mass conservation Eq. (1), while within the concentration boundary layer, the particle wall flux is determined with a one-dimensional semi-empirical particle deposition model (Lai and Nazaroff, 2000) and the results are substituted into Eq. (1) as the boundary condition.

Lagrangian approach

The equation of motion of a small aerosol particle can be written aswhere is the velocity of the particle, τ is the particle relaxation time, n(t) is the Brownian force per unit mass, ρ p and ρ are the particle and air density, respectively, and g is the gravitational acceleration. In a stochastically modeled turbulent flow, the instantaneous fluid velocity can be expressed aswhich is the sum of the mean velocity component, u, from the RNG k–ε model and the fluctuating velocity component, , which will be described in Section 2.3. The Brownian force per unit mass is important for submicron particles. The Brownian force is modeled as a Gaussian white noise random process as described by Li and Ahmadi (1992). The procedure for simulating the Brownian force is to generate a white noise process with the noise intensity,and the spectral intensity S 0 is given bywhere G is a zero-mean, unit variance-independent Gaussian random number and k B=1.38×10−23  J K−1 is the Boltzmann constant. The third terms on the RHS of Eq. (2) represents gravity force exerting on the particle.

It has been reported by the present authors that the requirement for a careful grid independence test is more stringent for a Lagrangian simulation (Lai and Chen, 2007). To model particle deposition accurately, the grid should be able to resolve not only the turbulence field, but also the flow field properties along particle paths. If the near-wall boundary layer is not properly resolved, the predicted deposition rate may depart from the actual solution severely. In the present model, the near-wall turbulence is deliberately damped. If the near-wall cell center is out of the viscous sublayer, the normal velocity at the point may not be negligible and consequently the interpolated normal velocity at the particle position may be sufficiently large to drive the particle to impact onto the wall directly. Hence, the near-wall grid should be fine enough to resolve the important deposition boundary layer. Three different grid systems are tested and details of them are listed in Table 1 . Taken into account of both computational resources requirement and result accuracy, grid system 2 with 181,976 hexahedral cells is chosen for the present simulation (Lai and Chen, 2007).

Table 1

Comparison of particle deposition fractions predicted with different grids

Particle size (μm)	Grid 1 (%)	Grid 2 (%)	Grid 3 (%)
7	61.7	59.1	58.7
1	48	11.3	11.1

Correction to the near-wall turbulent kinetic energy

To take into account of anisotropic turbulence, a new correction scheme, which is essentially a hybrid combination of the methods of He and Ahmadi (1998) and Matida et al. (2004), is proposed. The quadric relation used by He and Ahmadi is adopted, written aswhere A=0.008 (Bernard and Wallace, 2002) is used by fitting the DNS results of Kim et al. (1987). In order to implement this correction to the isotropic k–ε turbulence models, the method of Matida et al. is used to simplify the system by forcing the streamwise and spanwise normal Reynolds stress components equal to the normal component, i.e.A new turbulent kinetic energy for particle tracking calculations can then be defined as

In this method, only the turbulent kinetic energy in the near-wall area needs to be corrected, and the resultant turbulent time scale will be updated accordingly. It should be emphasized that the optimized turbulent kinetic energy in Eq. (8) still remains isotropic. There are several salient features of this hybrid scheme. The method can be conveniently applied to various conditions, as the modification to the model is minimal. The prediction of the normal Reynolds stress component is remarkably improved, thus it can result in a more reasonable particle deposition rate. On the other hand, the turbulent intensities in the rest two directions are underestimated undesirably. This under-estimation is unlikely to incur significant error when the particle deposition rate is the parameter of concern, as particle deposition is not directly affected by the rest two components and, the mean flow velocities in these two directions overwhelms the fluctuation counterparts.

Results and discussion

Fig. 2 depicts the particle concentration evolution of 0.3, 1 and 7 μm particles from time 0 to 1800 s (all particle sizes are referred to aerodynamic particle diameters). A noticeable observation is that the submicron particle concentration can achieve fairly uniform in just 300 s, and in 30 min the particle concentration is virtually homogeneous for the entire zone. The ultrafine particles exhibit similar dispersion characteristic to that observed for 0.3 μm particles (results not shown). The dispersion rate decreases with the particle size as concentration gradient is distinctively observed for coarse particles. For 7 μm particles concentration homogeneousness cannot be achieved (even beyond 1800 s).

Fig. 2

Typical concentration evolution of 0.3, 1 and 7 μm particles at four elapsed times, with inlet velocity 0.225 m s−1: (a) 0.3 μm, (b) 1 μm and (c) 7 μm.

To represent the non-uniformity of a concentration field, the coefficient of variation, , of concentration is defined as (Mage and Ott, 1996)where is the volume averaged concentration, C(t) is the concentration at elapsed time t at each sample point and N is the number of sample points. In case of non-uniform mixing, a useful parameter, the mixing time, is used to quantify the time needed to reach a well-mixed state. It is defined as the time at which the coefficient of variation becomes <10% permanently since the release of pollutant. Fig. 3 presents the coefficients of variation of the three particle sizes. Inferring from the results, it can be shown that CV strongly depends on the particle size and airflow; small size and high airflow favor mixing. One thing should be noted is that the mixing characteristics of 0.3 and 1 μm are very close whereas very distinct difference can be observed for 1 and 7 μm particles. As the gravitational settling magnitude is a function of particle size with an exponent of 2, sedimentation influences the mixing and dispersion rates (Fig. 2) in a non-linear way; in the results presented, the distinct behavior of 7 μm particles is observed.

Fig. 3

Coefficients of variation of concentration field for various particle size groups.

For Lagrangian model, the results presented are very different. Since particles are treated as a discrete and not a continuum phase, concentration contour does not exist. Another factor further complicates the issue is the computer resources required for Lagrangian model. Due to the current limitation of computational resources, almost all Lagrangian models inject particles momentarily and track the individual particle based on the forces exerted. Because of the transient in nature and the validity of the well mixed (discussed above), these pose limitation for presenting the results obtained from Lagrangian approach. In the literature, the parameter decay rate loss coefficient is often used to characterize the particle removal rate in an enclosure. The decay rate coefficient is derived from the mass balance principle where well-mixed condition is assumed. This has a significant implication on the application of particle decay rate coefficient for quantifying the fate of particles in a non-well-mixed enclosures and it must be used cautiously (Bouilly et al., 2005). Hence for Lagrangian framework, another common approach, which counts the number of particles remaining in the domain versus time, is adapted (Béghein et al., 2005). The particle deposition fraction, η, is defined aswhere N in and N out are the numbers of particles released at the inlet and the number of particles exited from the outlet, respectively. As observed in Fig. 4 , the original EIM without correction over-predicts particle deposition rate noticeably, particularly for 0.3 and 1 μm particles at the high ventilation rate. The turbulent intensity is stronger for the high-velocity case. If the correction scheme is not applied, the original EIM tends to give more severe over-prediction for the high-velocity case. Generally, the deposition fractions predicted with the correction are closer to the semi-empirical equation (Lai and Nazaroff, 2000).

Fig. 4

Deposition fractions to different surfaces predicted by the Lagrangian simulation.

Thus far we have highlighted the main features of two models developed recently and results are presented separately. It is very valuable to quantify the simulation results by a same parameter. In Fig. 5 , Eulerian and Lagrangian results for two airflow velocities are shown. Overall speaking, the results modeled by the two approaches agree well with each other; as the particle size increases, the deposition fraction increases. For submicron particles, the deposition fractions predicted by Lagrangian (without near-wall turbulent correction) is higher than those predicted with correction and Eulerian drift flux prediction follows. Deposition fraction for 0.45 m s−1 is also higher than the case of 0.225 m s−1. This is because for submicron particles, the deposition is significantly affected by turbulent diffusion and it is stronger for the high-velocity case. For supermicron particles scenario, the discrepancy between Eulerian and Lagrangian increases with particle size, however, is not significant. One thing should be noted is that as the particle size increases, the influence of the turbulent correction diminishes. These observations can be explained by the dominant fate mechanism for this size range of particles. As particle size increases, deposition by gravitational settling increases significantly and the other effect, i.e. turbulent diffusion, becomes less important. It must be careful when direct comparison is attempted to be made. As discussed above that the injection type for Eulerian is continuous while for Lagrangian, the injection is momentarily. For steady-state condition such as Eulerian model, the parameter deposition velocity can be used to quantify the deposition rate in a more rigorous way, while for Lagrangian framework, using deposition velocity as a parameter is not an obvious means. On the other hand, the deposition fraction cannot be used to quantify the deposition onto various orientated surfaces, however, it can be used to characterize both models.

Fig. 5

Comparison of Eulerian and Lagrangian models for particle deposition fraction.

The current study focuses on an empty chamber so some comments on practical applications should be addressed. The author has published one article regarding the effects of room furnishing on particle deposition rates indoors (Thatcher et al., 2002). The experimental results showed that increasing the surface area from bare to fully furnished increased the deposition loss rate with the latest increase seen for 0.5 μm particles. It can be attributed to the additional surfaces available for diffusion loss. The current models can be applied to both empty and furnished room with proper boundary conditions.

Conclusions

Knowledge of particle dispersion and deposition indoors improves human exposure assessment and provides insights for better pollutant control measures. CFD has become a virtual tool for studying particle dynamics. We compared a new drift–flux Eulerian and a modified Lagrangian models for a single-zone chamber geometry. One key feature of the drift–flux model is to encapsulate external forces into the formulation and present the results in a continuous domain. Near-wall anisotropic turbulence correction is properly applied in the Lagrangian model. Inferring from the results, it is shown that non-corrected turbulent scheme over-predicts particle deposition significantly for submicron particles while supermicron particles are not sensitive to the correction scheme. The paper also highlighted the results obtained by Eulerian and Lagrangian models. The two models agree very well for submicron particles. As the size increases, the discrepancy increases moderately. It is not straightforward to conclude which model is better than the other as each model has its own fundamental assumptions imposed. The results shown here reveal that the two present models are comparable.