Considerations for Modeling Particle Entrainment into the Wake of a Circular Cylinder

The objective of this work is to evaluate the performance of the steady state Reynolds Averaged Navier-Stokes (RANS) computational fluid dynamics (CFD) models for estimating concentration of low Stokes number aerosols (Stk = O(10⁻⁴)) in the wake of a bluff body. These simulations are compared with experimental data. In the simulations and experiments, particles are released upstream of the body and convected downstream, where some are entrained into the wake. The air velocity is computed using a steady state renormalized group k ∼ ϵ model. Lagrangian particle trajectory simulations are performed in conjunction with each airflow model to calculate concentrations. The experiments are performed in an aerosol wind tunnel in which phase Doppler velocimetry measurements are obtained for the velocity field and aerosol concentration.

The RANS model yields a wake concentration deficit that extends downstream past x/D = 10, while the experiments produce elevated concentrations immediately downstream of the near wake. It is postulated that the concentration peak is at least in part attributed to particle interaction with the boundary layer by the following mechanism. Particles are transported into the boundary layer by turbulent diffusion, turbophoresis, and/or inertial forces. Particles then separate from the cylinder with the airflow and travel in a sheath around the periphery of the near wake to converge at the downstream edge of the near wake. Underestimation of the wake concentration by the RANS model is potentially due to inadequacy in the boundary layer approximation used in the model.

INTRODUCTION

The preponderance of human exposure computational fluid dynamics (CFD) studies employ Reynolds Averaged Navier-Stokes methods (RANS) (e.g., Gadgil et al. 2003; Welling et al. 2001; Lan and Viswanathan 2001; Kulmala 2000; Flynn and Sills 2000). This choice may be related to the availability of RANS models in commercial CFD packages. Many software packages allow easy creation of computational meshes to accommodate more complex geometries. Furthermore, the coefficients used in RANS models are generally considered robust enough to fulfill a wide range of engineering problems (Launder and Sharma 1974; Bradshaw 1997). It is also much less costly to run a RANS simulation than a large eddy simulation or a direct numerical simulation (Bradshaw 1997). For the engineer, who often requires an answer in terms of average phenomena, RANS would seem to be a sensible choice.

For bluff body flows, shortcomings of k ∼ ϵ models have been well documented. These have been attributed to a variety of reasons, including:

•	inability of these fluid models to capture boundary layer behavior (Majumdar and Rodi 1985);
•	inaccuracies in the boundary layer approximation for turbulence kinetic energy (Matida et al. 2000);
•	omission of periodic vortex shedding (Lübcke et al. 2001; Majumdar and Rodi 1985; Gomes et al. 1997);
•	failure to represent the interaction between Reynolds stresses and mean flow (Lübcke et al. 2001);
•	late boundary layer separation (Rodi 1997; Bosch and Rodi 1996; Bosch and Rodi 1998; Kim and Boysan 1999);
•	excessive turbulence kinetic energy at the upstream stagnation point of the bluff body, leading to underestimation of the near wake size (Chen and Kim 1987; Jenne and Reuss 1999; Henkes 1998);
•	overprediction of turbulence dissipation (Shih et al. 1995; Lumley 1992; Yakhot and Orszag 1986); and,
•	inaccurate turbulence closure approximations (Franke et al. 1989; Mompean et al. 1996).

These factors have resulted in poor validation with experimental data and results from large eddy simulations (Richmond-Bryant 2003a; Lübcke et al. 2001; Breuer 1998). It can be postulated that use of RANS models in conjunction with a particle simulation would necessarily cause inaccuracies in aerosol trajectories and concentrations. However, this requires further proof. The work presented here uses the RANS model not because it provides the most accurate representation of the near-wake physics. It is applied because it is so commonly used in applied fluid dynamics simulations that the effect of this model on a near-wake aerosol concentration field should be demonstrated in the literature.

The most common approach to modeling aerosol trajectories in CFD is to use an Eddy Interaction Model (EIM) (Gosman and Ioannides 1981). This is given by the following group of equations:

where τ = aerosol relaxation time, u _p = aerosol velocity, u = air velocity, and F _p = external forces on the aerosols. The air velocity is represented in EIM as the sum of mean and fluctuating velocity components:

where λ = random number generated from a uniform distribution. In EIM, λ remains constant during one lifetime of the large-scale eddy that the aerosol is crossing, although it has been suggested that two eddy lifetimes is more appropriate for preserving proper dispersion (Graham and James 1996). Given the dependence of Equation (1) on Equation (2), the turbulence kinetic energy becomes very influential on accuracy of the particle trajectory.

Studies in the literature have demonstrated that use of the EIM in conjunction with boundary layer approximations to be deficient for small particle transport in pipes (Matida et al. 2000). This is largely due to the contribution of the turbulent velocity fluctuations on the overall fluid velocity in the aerosol transport equation. This has been addressed for aerosol transport in the boundary layers of pipes using a semi-empirical formulation for the turbulence kinetic energy in the boundary layer (Matida et al. 2004). This is in lieu of the normal boundary layer model approach of assuming that the turbulence kinetic energy at the freestream edge of the boundary layer cell is constant throughout the cell (Launder and Spalding 1974). In reality, turbulence energy approaches zero with decreasing distance from a boundary. Given that the boundary layer model is already based on an empirical approach (Reichardt 1951), it is sensible to enhance that model with an empirical formulation for the turbulent fluctuations. However, it is unclear that the boundary layer turbulence model used by Matida et al. (2004) is appropriate for anything other than turbulent pipe flow because the structure of these flows is quite different. Furthermore, given the frequency of applying RANS and EIM in conjunction with bluff body studies of aerosol transport (e.g., Gomes et al. 1997; Hu and Hsiao 2005), it is appropriate to examine the results from a RANS/EIM simulation in light of experimental data. Therefore, this study explores the effect of EIM with the wall model on bluff body flows prior to implementing a modification to the methodology.

In the present paper, we examine the impact of the wall model on the transport of fine particles (particle diameter, d _p , = 3.5 μ m, particle density, ρ _p , = 2400 kg/m³) near bluff bodies. The time-averaged concentration profile of a dilute suspension of fine particles within the wake of a two-dimensional circular cylinder obtained with a RANS fluid model and EIM particle tracking model is compared with experimental data. Particles are released upstream of the cylinder to permit observation of average particle entrainment into the recirculation zone downstream of the cylinder. Particle concentration is computed at discrete time steps from the velocity field and time-averaged.

NUMERICAL METHODS

RANS Airflow Model

The two-dimensional, turbulent Navier-Stokes equations are solved using the renormalized group (RNG) k ∼ ϵ turbulence model with a Boussinesq constitutive relationship. The system of equations is solved using FIDAP v.8.6.2 on an IBM 933 MHz processor PC, an SGI Origin 2400 with 48–400 MHz processors, and a Sun workstation with Ultra Sparc III dual 900 MHz processors. Details of the airflow simulation are described in Richmond-Bryant (2003a) and are summarized here. To take advantage of steady-state symmetry, the velocity field is solved over a half-plane, semi-circular grid. Boundary conditions are selected to match wind tunnel experiments of the airflow (Richmond-Bryant 2003b). All parameters in the airflow simulations are nondimensionalized. Laminar viscosity, μ, is set equal to 1/Re, where Re = Reynolds number, Re = ρ DU _∞ /μ = 5,232, U _∞ = freestream velocity, and D = cylinder diameter. Inlet turbulence kinetic energy is calculated using the measured freestream turbulence intensity, TI, in the tunnel: k = 1.5 (TI· U _∞)² (Wilcox 1994). With TI = 3.5%, the freestream turbulence kinetic energy is 1.84 × 10^{− 3}. From the relationship between turbulent viscosity, turbulence kinetic energy, and dissipation, a relationship for ϵ is derived: ϵ = C _μ v ρ k ²/R _μ μ, where R _μ = μ _t /μ(Launder and Spalding 1974; Bosch and Rodi 1998). Based on Bosch and Rodi (1998), it is assumed that R _μ = 10 for the ϵboundary condition; this yields an inlet boundary condition of ϵ = 1.59 × 10^{− 4}. Outlet boundary conditions are left free. At the cylinder, no-slip conditions apply and = 0. Along the symmetry plane, there is zero tangential stress and normal velocity. Verification testing of the simulation follows Roache (1998) and is described in detail in Richmond-Bryant (2003a).

Figure 1 shows a comparison between the experimental and simulated streamline velocity field for three x/D positions downstream of the cylinder. The simulations clearly capture the shape of the velocity field in the wake. However, these results demonstrate that the simulations underestimate the width and length of the wake. This is an expected result, considering the findings of (Bosch and Rodi 1996, Bosch and Rodi 1998), Kim and Boysan (1999), and others. As mentioned in the introduction, the observed discrepancy raises the question whether it is valuable to pursue a comparison of the experimental and computational aerosol concentration fields in the near wake. However, given the paucity of experimental and numerical comparisons of near-wake concentration fields for the still commonly used k ∼ ϵ models, this demonstration is worthwhile for the applied modeling audience.

FIG. 1 Cross-stream velocity field at various axial distances from the cylinder for the experiments and simulations. The cylinder is centered at (0,0).

RANS Aerosol Transport Modeling

Particle transport in FIDAP is based on the EIM proposed by Gosman and Ioannides (1981) for inserting into the particle drag equation the point-wise air velocity, u, which is a function of mean and fluctuating components. New developments in the EIM for using an empirical formulation for the turbulent velocity fluctuations in the boundary layer were not used here because they are based on data for turbulent pipe flow (Matida et al. 2004). It has not yet been shown if this formulation is appropriate for bluff body transport of particles.

For these simulations, the error resulting from the time-step size is assessed by repeatedly releasing one particle in laminar flow with subsequently smaller time-steps. Laminar flow is designated to evaluate the error without the variability resulting from turbulence. The time-step is halved until Δ t = 0.003125, when errors in x _p , y _p , u _p , and v _p are reduced to O(10⁻³). In the interest of computational time, final simulations are performed with Δ t = 0.00625 since the errors between results for the Δ t = 0.00625 simulation results and the Δ t = 0.003125 simulation results are very small.

It is assumed that the cylinder has a perfect bounce boundary condition to emulate the experiments using a smooth aluminum cylinder and solid particles. This is modeled in FIDAP with a restitution coefficient of 1 on the cylinder. In reality, there is probably some deposition of the particles onto the cylinder to yield a restitution coefficient slightly less than 1. This assumption is based on Dahneke's (1971) analysis of aerosol impaction with flat and cylindrical surfaces, in which he noted that any plastic deformation occurring during particle-surface interaction is improbable for aerosol particles because the size of the particle is likely on the same scale as small surface irregularities. Likewise, Dahneke's (1971) calculations show that presence of a cylinder in place of a flat plate should be inconsequential when the cylinder diameter is orders of magnitude larger than the particle, as is the case for our study.

Final particle simulations are performed with the turbulence model invoked. 3.5 μ m particles with ρ = 2,400 kg/m³ are released when t = 0 at 101 evenly spaced points along the upstream line segment 0 ≤ y/D ≤ 1.5 at x/D = −5. These quantities are nondimensionalized for the simulation through Stokes number conversion. The resulting trajectories are added to a dataset generated from all previous runs and input into a particle concentration algorithm, described below.

RANS Particle Concentration Field

To compute particle concentration, a rectangular segment of the domain from (−5, 0) to (10,3) is discretized into 3025 square cells with dimension 0.125 × 0.125. During each time-step, particle concentration is computed in each cell. Particle number concentration is adapted from Heinsohn's (1991) formula for particle mass transport through a defined region of air:

where A _c = area of one cell (equal for all cells), t _i = length of time the particle spends within the cell during the time-step, n _gi = number generation rate for the ith particle trajectory, and N = total number of particles crossing the cell. Particle concentration is computed using an in-house code taking input from particle trajectory files generated by FIDAP. It should be emphasized that the cells in which concentration was computed do not correspond to the cells used for finite element solution of Navier-Stokes. The concentration cells are designated after solution of the velocity field in a post-processing step. Concentration values for each time-step are factored into the time-averaged concentration for each concentration cell. Concentration cell length scale, L, is chosen based on Courant criteria, adjusted for turbulent velocity fluctuations, as L = 0.125 ≫ (U _max+ λ U′_max)Δ t ≈ 0.0306. 202,000 particles are used in the simulation. This number of particles has been selected based on Graham and Moyeed's (2002) statistical strategy for determining if the number of particles in a Lagrangian simulation is sufficient for convergence of the concentration field.

EXPERIMENTAL DATA

Experiments were performed in a recirculating wind tunnel located at the U.S. Environmental Protection Agency National Exposure Research Laboratory (Research Triangle Park, NC). The tunnel consisted of a blower powered by a 50-hp motor that drew air downstream through a set of louvers. Prior to entering the test section, the air flowed through a honeycomb grid, which removed any large-scale turbulent structures from the airflow. Then, air traveled downstream to the 1.5 m × 1.2 m × 7.3 m test section, where velocity measurements were taken with a two-dimensional phase Doppler anemometer (PDA) (Dantec Measurement Technology, Inc., Copenhagen, Denmark). For a freestream velocity of 1.0 m/s, the tunnel velocity profile was reasonably uniform over the cross-section of the test area, with a spatial coefficient of variation (CV) of 3.0%. The turbulence intensity at this speed was measured at 3.5%.

Particles were aerosolized by a fixed array of Venturi nozzles fed by a conveyor belt coated with particles. A detailed description of this aerosol generation system can be found in Heist et al. (2003). The aerosol generation system yielded a spatially uniform concentration profile across the wind tunnel test section with a CV = 6.4%. Nephelometry was used to demonstrate that the particle concentration remained relatively constant in time with a CV = 6.2% over the 2-minute sampling period.

In these experiments, a polydisperse aluminum silicate dust (Zeeospheres^TM W-610, 3M, St. Paul, MN) with a count mean aerodynamic diameter of 5.16 μm (count mean physical diameter = 3.33 μm) and density of 2400 kg/m³ was released upstream of the honeycomb. The aerosol was allowed to mix with the air before traveling downstream past an infinite 0.0762 m-diameter cylinder mounted in the test section. To minimize charging effects, the cylinder was grounded for all experiments. PDA measurements of particle velocities were taken at 102 locations in the cylinder wake and converted to point-wise concentrations for each particle size. Because they were most prevalent in each sample, 3.5 μm particles (physical diameter) were used in the data analysis. Additionally, to measure the air velocity field, the air was seeded with propylene glycol fog particles (Martin Magnum Pro 2000, Århus, Denmark) in separate experiments and measurements were taken at 140 locations in the cylinder wake (Richmond-Bryant 2003a).

The phase Doppler anemometer used for testing the air and particle velocity profiles had a multiline Argon-Ion laser with 488 nm and 514.5 nm TEM₀₀ beams (Coherent, Model Innova 300, Santa Clara, CA) for the U and V components of velocity, respectively. These were passed through a transmitter (Dantec, Model 60X41, Skovlunde, Denmark) containing a Bragg cell to create a frequency shift of 40 MHz in each of the two beams. The beams were focused with a separation width of 480 mm for the U velocity component and 472 mm for the V velocity component, and the separation angle was 22°. This configuration yielded a fringe spacing of 1.4576 μ m for the U velocity component and 1.5119 μ m for the V velocity component. There were 32 fringes in the x-direction and 33 in the y-direction. The signal produced by refraction of light when particles pass through the measuring volume was captured by photomultiplier tubes (Dantec, Model 57X08, Skovlunde, Denmark), registered by the signal processor (Dantec, Model 57N10, Skovlunde, Denmark), and then sent to a computer for statistical analysis of the velocity components, particle size distribution, concentration levels, and turbulence intensity.

Prior to analyzing the data, provisions for minimizing biases and uncertainties related to velocity measurements and concentration computations were made. Instrumentation bias resulting from trajectory ambiguity and slit effect (Maeda et al. 1996) were minimized through hardware adjustments of Bragg cell frequency shifting and the use of a third detector to obtain two phase differences by which the particle can be sized (Maeda et al. 1996; Bates and Ayob 1995). Uncertainties in the particle velocity field related to signal-to-noise ratio and sample independence were relatively low, O(10^{− 3}). These velocity uncertainties were also used in creation of confidence intervals for the experimental concentration. Details of each computation are given in Høst-Madsen and Caspersen (1994) and George (1979) for signal-to-noise ratio and Winter et al. (1991) for sample independence. Based on the optical constraints, limit of detection should not be an issue because laser Doppler anemometry is a primary standard measurement. This is based on the contingency that the photodetectors are aligned for the optimal scattering angle of the seed particles. However, data are subject to machine round-off error in signal processing. For the 8-bit signal processor, the round-off error results in computation to three significant digits and the limit of detection is 0.001 m/s.

Particle concentration is computed under the assumption that each particle contributes to the length of the measurement volume (Qiu and Sommerfeld 1992; Sommerfeld and Qiu 1995; Aísa et al. 2002):

where T = total sampling time at a point, A = effective area of the measurement volume, U = measured velocities of each individual particle, N = total number of particles. The subscript j refers to individual particles, and the subscript k represents particles moving in the direction. In this analysis, only the x-direction is used for k. Equation (4) is based on the assumption that the particles move in one predominant flow direction but accounts for the small number of particles that may deviate from the freestream. Hence, it is more applicable for turbulent flows than are formulations using the average velocity of the particles. Because there are limitations in the PDA equipment available for this study, more current signal processing methods to account for directionality, such as vector volume flux (Berkner et al. 2001) or integral methods (Sommerfeld and Qui 1998) are not employed here.

The primary source of uncertainty for computing particle concentration in the near wake of a cylinder is determining whether all counted particles are distinct or if one particle travels back and forth through the sampling volume to generate multiple counts. If all counted particles are distinct, then Equation (4) is directly applicable. If the same particle travels multiple times through the sampling volume, then the flux of the particles moving against the predominant direction of transport is subtracted from the flux of those moving with the predominant direction. When most of the particles move in one direction, these two approaches yield close results, and the uncertainty is relatively low. When roughly half the counted particles move in one direction and half move in the opposite direction, then the difference between these two quantities increases and the uncertainty is high. The two different methods can be used to estimate upper and lower bounds for the concentration profile. Using this methodology, regions of high uncertainty are identified to exist in the near wake (x/D < 2.5, 0 ≤,y/D ≤ 0.5). Downstream (x/D ≥ 2.5) and along the sides of the near wake (y/D > 0.5), particle transport is predominantly unidirectional and this uncertainty is reduced. This is consistent with Saffman's (1987) observations that directional bias must be employed to correct the effective area of the measurement volume. Saffman (1987) uses the average angle of approach for this purpose. A large directional bias correction term generally means that there is a lot of variation in angle of approach for each particle. For this reason, sampling points where the adjusted area is greater than 110% of the area computed by the software are discarded from the comparison with simulation results.

RESULTS

Figure 2 compares the time-averaged cross-stream aerosol concentration field obtained with the steady RANS simulation with that obtained through the experiments averaged over a 2-minute period (43.5 vortex shedding cycles). Clearly, the RANS simulation underestimates the concentration in the near wake and at the downstream edge of the near wake (around x/D = 2.5), even with the large experimental uncertainty bounds within the near wake. These concentration fields are shown graphically in Figure 3 and Figure 4 for the experiments and simulations, respectively.

FIG. 2 Cross-stream aerosol concentration field at various axial distances from the cylinder for the experiments and simulations. Concentrations are normalized by the freestream aerosol concentration upstream of the cylinder, and the cylinder is centered at (0,0).

FIG. 3 Experimental concentration profile map for the domain around the circular cylinder, which is centered at (0,0) comparing a) the lower bound (top) and b) the upper bound (bottom).

FIG. 4 Simulated concentration profile map for various distances downstream of the cylinder, which is centered at (0,0). All concentrations are normalized by an upstream freestream concentration so that a concentration of 1 is equivalent to that of the freestream.

Upper and lower uncertainty bounds on the experimental measurements of concentration are computed and displayed with a map in Figure 3. These intervals are based on the estimated uncertainty from sampling in turbulence, signal noise, and sample independence; they do not represent 95% confidence intervals. There are large discrepancies between the lower and upper concentration bounds within the near wake region. This highlights the high level of uncertainty in the near wake. The upper bound shows a large surplus in the near wake concentration, while the lower bound of concentration appears to be less than one at many locations. Along the outer edge of the wake (y/D > 0.5) and downstream of the near wake (x/D > 2.5), particle transport is more nearly unidirectional and the uncertainty decreases drastically. The concentration profiles in the upper and lower bounds are nearly identical in these regions of the domain. A peak in concentration is seen immediately downstream of the near wake at 3 ≤ x/D ≤ 3.5 over the span ∣ y/D∣ ≤ 0.5. The maximum in this peak occurs around x/D = 3, y/D = 0.375. Downstream of this peak, the concentration returns to freestream levels. Significant elevated concentrations are also seen around the wake edge at x/D = 2, y/D = 0.375.

For the steady-state concentration profile obtained using the RANS simulation shown in Figure 4, a concentration deficit exists for the near and far wakes. On the upstream side of the cylinder, a high concentration sheath can be seen bordering the cylinder from the stagnation point to approximately 80° around the cylinder. At this point, concentration levels bordering the cylinder return to freestream levels prior to the separation point, computed at 131°. Mass balance across various x/D planes demonstrates that 95% of trajectories are conserved at x/D = 5. Occasional losses occur at the symmetry plane, by retention near the cylinder, and in the near wake.

DISCUSSION

Identification of the elevated concentrations around the wake and near the average downstream stagnation region obtained in the experiments is not intuitive. However, these high concentration regions have also been observed by Fohanno and Martinuzzi (2004) and Jacober and Matteson (1990), who also demonstrate C/C₀ > 1 for small particles at the edge of the wake and downstream of a bluff body. Elevated particle concentrations near boundaries have also been shown experimentally by Sun and Lin (1986) and Sehmel (1968) and numerically by Matida et al. (2000) for small particles transported in turbulent pipe flow.

When considering the possible effect of boundary layer dynamics on the particle concentration field, it is useful to keep in mind the inertial characteristics of the particles in our investigation when considering developments in the literature. Based on the particle size characteristics in our study (d _p = 3.5 μ m, ρ _p = 2400 kg/m³) it is estimated that T ⁺ _max∼ 0.2 in the boundary layer. Here T ⁺ = the nondimensional particle relaxation time, T ⁺ = T _p / T _e . In the boundary layer, the eddy timescale is developed from the friction velocity and the kinematic viscosity, v: T _e = v/u*². The estimated eddy timescale for computing our T ⁺ is obtained using data on friction velocity from the shear stress profile around a circular cylinder at Re = 3,900 computed with a large eddy simulation (Breuer 1998). Although T ⁺ = 0.2 is not particularly large, it is sufficient for particles to exhibit some inertial behavior. If the high concentration region near the stagnation point shown in the simulation is true of the experimental results, inertial deviation from the fluid streamlines may play a minor role.

Young and Leeming (1997) analyzed Reynolds averaged equations of conservation of particle mass and momentum when studying particle deposition in pipes. In doing so, they also gained some perspective on the particle concentration near the pipe wall. In this work, they found that two important effects influence particle deposition near walls: diffusion and turbophoresis. Based on this and previous studies of particle deposition in channels, they defined three regimes in which the particles can be classified with respect to how they are influenced by turbulence. When T ⁺ < 0.3, particles fall in the diffusional deposition regime. In this regime, Young and Leeming (1997) found that primarily “particles are transported into the sub-layer by turbulent diffusion and finally to the wall by Brownian diffusion” given the correlated fluctuations of the particle location in the x- and y-direction under turbulent diffusion since u'v′ > 0 where > 0 over the boundary layer (Hinze 1975). Hence, turbulent diffusion is chiefly responsible for the build-up of particle concentration around the boundary layer. Particles then could deposit on the walls via Brownian diffusion.

Turbophoresis was demonstrated by Young and Leeming (1997) to have some impact on particle concentration for T ⁺ = 0.2, although this effect was thought to be much less significant than turbulent diffusion. Our estimation of the turbophoretic velocity of a particle, given nondimensionally as (Reeks 1983):

with T ⁺ = 0.2 shows that it is, at most, ∼ 1% of the mean air velocity in the boundary layer in regions where the wall model is applied. Here, ⟨ k ⟩ is the averaged nondimensional turbulence kinetic energy. Near the upstream stagnation point, however, the wall model is not valid. In this region, turbophoretic velocity becomes significant in comparison with the near-zero mean boundary layer velocity. In this manner, particles collecting near the stagnation point via turbophoresis could then diffuse into the boundary layer. Particle retention at the boundary layer edge may still be controlled by turbulent diffusion in locations along the cylinder wall between the stagnation and separation points.

Figure 5 displays a close-up of the RANS simulated concentration profile around the cylinder boundary. Recall that a concentration of 1 is equivalent to the upstream freestream concentration after normalization. A sharp peak in concentration occurs at the upstream stagnation point on the cylinder. This is consistent with our hypotheses of turbophoretic induction of particles at a location of zero mean velocity and inertial deviation of particles into the stagnation point. The high concentration sheath around the cylinder is maintained on the upstream side up to θ ∼ 80°. At this point, the sheath disappears. It is hypothesized earlier in this discussion that particles remain in the boundary layer under the influence of turbulent diffusion until the separation point, modeled at 131° by the FIDAP airflow simulation. This presents the question of how the model handles turbulent diffusion in the boundary layer and what occurs around θ = 80° that could possibly prevent this mechanism. Although the FIDAP model yielded late separation, as shown in Figure 1 and in the literature (Rodi 1997; Bosch and Rodi 1996; Bosch and Rodi 1998; Kim and Boysan 1999), this error should only affect the size of the wake region and not cause a deficit of particles in the far wake. The unphysical assumption of constant turbulence kinetic energy in the boundary layer causes overestimation of velocity fluctuations responsible for turbulent diffusion in the boundary layer. As already indicated in Matida et al. (2000), this could lead to particles being ejected from the boundary layer and, hence, underestimation of the particle concentration in the boundary layer.

FIG. 5 Close-up of particle concentration profile around the cylinder surface. The upstream freestream concentration is used to normalize all concentrations so that a concentration of 1 is equivalent to that of the freestream.

To examine this potential source of error, the turbulence kinetic energy profile around the edge of the cylinder is examined, as shown in Figure 6. Around θ = 80° a sharp increase in turbulence kinetic energy occurs. This location roughly coincides with the disappearance of the high concentration sheath. It is possible that the sharp increase in turbulence kinetic energy causes an overestimation of particle fluctuation in the near-wall region. This would then lead to the ejection of particles from the near-wall region around θ = 80°. Clearly, use of an empirical formulation for the turbulence kinetic energy around the cylinder edge like that used in Matida et al. (2000) for pipe flows would yield improvements that may prevent the particles from being ejected from the boundary edge. However, differences in the boundary layer velocity and turbulence kinetic energy profiles in regions of upstream stagnation, separation, and in between would present an argument for detailed measurement of the boundary layer velocity and energy profiles around a circular cylinder for use of empirical boundary layer formulations with CFD modeling.

FIG. 6 Turbulence kinetic energy profile around the cylinder boundary. Note that the turbulence kinetic energy is constant over the height of the cells bordering the cylinder boundary.

CONCLUSIONS

In summary, comparison between experimental results measuring aerosol concentration in the wake of a circular cylinder and results from RANS simulations demonstrate large discrepancies in the resulting time-averaged concentration profile for 3 μm aerosols released upstream of a cylinder. While the steady RANS model produced a concentration deficit throughout the wake, the experimental results demonstrated a concentration surplus around the wake edge and further downstream. We explore the possibility that use of a steady simulation and application of the wall model in FIDAP software cause these differences. Inaccurate modeling of the turbulence kinetic energy profile in the boundary layer may produce overestimation of turbulent particle diffusion around the cylinder wall. These errors may then promulgate to low predictions in the wake concentration profile. Future work should follow Matida et al. (2000) with experimental development of the turbulence kinetic energy field in the boundary layer around a circular cylinder for use in an empirical boundary layer model in CFD simulations.

Acknowledgments

The research presented in this article was supported in part by grants #1 R01 OH07363 and #1 R03 OH09302-01 from the National Institute for Occupational Safety and Health (NIOSH), Centers for Disease Control (CDC). Its contents are solely the responsibility of the authors and do not necessarily represent the official views of NIOSH. Computational resources were made available by the North Carolina Supercomputing Center, project 11059.