Skip Navigation


Annals of Occupational Hygiene Advance Access originally published online on September 12, 2005
Annals of Occupational Hygiene 2006 50(2):157-173; doi:10.1093/annhyg/mei048
This Article
Right arrow Abstract Freely available
Right arrow FREE Full Text (PDF) Freely available
Right arrow All Versions of this Article:
50/2/157    most recent
mei048v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Add to My Personal Archive
Right arrow Download to citation manager
Right arrowRequest Permissions
Google Scholar
Right arrow Articles by ANTHONY, T. R.
Right arrow Articles by FLYNN, M. R.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by ANTHONY, T. R.
Right arrow Articles by FLYNN, M. R.
Social Bookmarking
 Add to CiteULike   Add to Connotea   Add to Del.icio.us  
What's this?


© 2005 British Occupational Hygiene Society Published by Oxford University Press


Original Article

CFD Model for a 3-D Inhaling Mannequin: Verification and Validation

T. RENEE ANTHONY1,* and MICHAEL R. FLYNN2

1 Environmental and Community Health, Mel and Enid Zuckerman College of Public Health, The University of Arizona, Tucson, AZ 85721-0468, USA; 2 Department of Environmental Sciences and Engineering, School of Public Health, The University of North Carolina, Chapel Hill, NC 27599-7431, USA

* Author to whom correspondence should be addressed. Tel: +1 520 8825852; fax: +1 520 882 5014; e-mail: tra{at}email.arizona.edu


    ABSTRACT
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
This work investigates the use of computational fluid dynamics (CFD) to model air flow and particle transport associated with an inhaling anatomical mannequin. The studied condition is typically representative of occupational velocities (Re = 1920) and at-rest breathing (R = Uo/Um = 0.11). Methods to verify and validate CFD simulations are detailed to demonstrate convergence and describe the model's uncertainties. The standard k-epsilon model provided a reasonable flow field, although vertical velocity components were consistently smaller than the experimental validation data, owing to truncation of the computational model at hip height. Laminar particle trajectory studies indicated that the modeled velocity field resulted in a shift of particle aspiration fractions toward particles smaller than those determined experimentally, consistent with the vertical velocity field differences.

Keywords: inhalable particles • modeling computational fluid dynamics


    INTRODUCTION
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
Particle inhalability, or aspiration efficiency, is defined by the ratio of the inhaled particle concentration to the uniform freestream particle concentration. This ratio varies with particle size and is affected by the freestream velocity, breathing rate and pattern and orientation of the breather relative to the flow field. Traditional studies to investigate these ratios rely on wind tunnel experiments where mannequins inhale particles uniformly distributed in the wind tunnel. However, as particle sizes increase and freestream velocities decrease, the ability to ensure spatially uniform particle distributions becomes difficult. Hence, the denominator in the particle inhalability calculation becomes more uncertain. Computational fluid dynamics (CFD) modeling of particle aspiration circumvents this fundamental problem inherent in wind tunnel experiments.

Numerical simulations have been used for the past 24 years to explore particle inhalability (Table 1). Initial studies used simple, 2-D shapes to represent the human form (Ingham, 1981Go; Dunnet and Ingham, 1986Go; Ingham and Hildyard, 1991Go; Chung and Dunn-Rankin, 1992Go). With advances in computer speed and memory, researchers performed 3-D simulations by representing the human head as a sphere (Dunnett and Ingham, 1988Go; Dunnett, 1997Go, 1999Go; Dunnett and Vincent, 2000Go) with a circular opening for the mouth. However, these models ignore the bluff body effect of the human torso, in which the velocity of the approaching freestream under the head decreases streamwise but increases laterally and vertically. Erdal and Esmen (1995)Go added a torso effect by using a circular cylinder with a hemispherical top to investigate particle aspiration through a circular mouth opening using potential flow simulations. To date, these models typically examined high-velocity conditions rather than the low-velocity (<0.3 m s–1) freestream representative of occupational environments (Baldwin and Maynard, 1998Go).


View this table:
[in this window]
[in a new window]
 
Table 1. Review of existing computational simulations of particle inhalability

 
Dunnet and Ingham (1988)Go used the boundary integral equation method to investigate aspiration efficiency of a human head by studying inhalation of a sphere dimensioned to match the mannequin head used in the wind tunnel experiments of Ogden and Birkett (1977)Go. They found that the sphere provided similar results to the head studies, except at low suction (suction velocity = 1.07 m s–1, freestream velocity = 2.75 m s–1), where the facial features created greater effects on the fluid flow near the mannequin. Recent experimental work comparing aspiration differences between simple geometric and more complex anatomical mannequins indicated that, in low freestream environments, both the velocity fields and the concentrations of aspirated particles were significantly different between these mannequins (Anthony et al., 2005Go). The authors recommended using a 3-D model with facial details to investigate particle inhalability in CFD studies.

Other recent CFD investigations of flow past heated mannequins have successfully used more realistic shapes to represent humans (e.g. Hyun and Kleinstreuer, 2001Go; Hayashi et al., 2002Go), an improvement over the cubic-mannequin forms and estimated inhalation volume used in Brohus and Neilsen (1996)Go and Brohus (1997)Go. However, these studies are of limited use to particle inhalability studies owing to the unrealistic breathing simulations. Murakami et al. (1999)Go ignored breathing and used a Launder-Sharma low-Reynolds k-epsilon model to examine the thickness of the boundary layer as a function of freestream velocity (0.25 and 2.5 m s–1) and turbulence intensity (6 and 29%). Hayashi et al. (2002)Go used a standard k-epsilon turbulence model with inhalation to examine the effects of body positions on the rising thermal plume. Because this model focused on room ventilation, it is of limited use to inhalability studies requiring uniform particle concentration and freestream conditions. Hyun and Kleinstreuer (2001)Go used the RNG k-epsilon turbulent model to study inhaled concentrations of tracer gases by the most realistic of humanoid shapes in these studies. Although they considered respiration, thermal heating, and room ventilation conditions, their mannequin head was devoid of facial features that affect particle transport and inhalation.

The goal of this work is to develop a CFD model to examine particle inhalability associated with a human-shaped mannequin, with emphasis on quantifying the numerical uncertainties. This paper investigates the fluid dynamics associated with an inhaling humanoid mannequin and discusses the adequacy of the standard k-epsilon turbulence model in low-velocity environments. Model verification includes a detailed discussion of the two major sources of computational error: non-linear iteration and mesh convergence (Roache, 1998Go; Stern et al., 2000Go). The uncertainties of numerical velocity prediction are assessed using the grid convergence index (GCI) (Roache, 1998Go). Experimental velocity measures are used to validate the numerical model by examining where the uncertainty bands of predicted and measured velocities overlap. Subsequently, both laminar and turbulent particle trajectories are studied to explore the particle transport model.


    METHODS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
The 3-D air velocity field associated with an inhaling mannequin facing the wind was measured using laser Doppler anemometry (LDA) in a wind tunnel and was subsequently simulated numerically using Fidap v. 8.7.2 (Fluent, Inc., Lebanon, NH). Details of the experimental conditions have been discussed (Anthony et al., 2005Go), but a summary is provided here for clarity. Figure 1 illustrates the geometry, coordinate system and nomenclature used throughout this study.


Figure 1
View larger version (30K):
[in this window]
[in a new window]
 
Fig. 1. Geometry, coordinate system, and nomenclature used throughout this study.

 
Experimental
The test mannequin, modified to match the anthropometric dimensions of a full-scale female, was a commercially available doll, My Size BarbieTM (Mattel Inc., El Segundo, CA). The mouth opening was shaped within the lips as a round-edged rectangle (0.01645 m wide by 0.00396 m tall, area = 3.09 x 10–5 m2). The origin of the coordinate system was at the mouth center for both numerical and experimental work. The mannequin was positioned in the wind tunnel to face the wind. The ratio of inhalation velocity (2.7 m s–1) to freestream velocity (0.3 m s–1) was 0.11. The Reynold's number, using the head diameter (0.096 m) as the length scale, was ~1910 and represents a person in a low velocity (0.2 m s–1) environment.

The investigation of the fluid flow field was conducted in a closed-loop wind tunnel, measuring 1.2 m in height and 1.5 m in width, with a 7 m glass-walled test section. The mannequin blocked <6% of the wind tunnel cross section. Velocity measurements were obtained using a LDA (Dantec 41N10 Traversing Amplifier, 58N10 PDA Signal Processor and Sizeware 2.0 software; Dantec Dynamics A/S, Solvlunde, Denmark). A theatrical smoke generator (Martin Magnum Pro 2000, Martin Professional, Denmark) was used to seed the flow for the measurements.

The LDA used a Coherent Innova 306 argon-ion laser (Laser Innovations, Moorpark, CA) to generate laser light that was split into the three wavelengths for velocity measures: 488.0 nm (blue, for streamwise velocity), 514.5 nm (green, for lateral velocity) and 476.5 nm (violet, for vertical velocity). The Dantec FiberFlow optical system split each of these beams in two, passing one of each pair through a Bragg cell, thereby shifting the frequency by 40 MHz to eliminate the ambiguity in velocity direction measurements. Doppler burst signals from the photomultiplier tubes were analyzed to obtain velocity measurements.

The 3-D, time-averaged velocity components were measured upstream of the inhaling mannequin. Measurements were taken in 0.005 m increments laterally and vertically to distances ±0.030 m and at distances of 0.011, 0.015, 0.020, 0.030, 0.040, 0.050 and 0.100 m upstream of the mouth opening for a total of 1171 locations. Because only half the velocity field was investigated in the numerical simulation, velocities on the two sides of the symmetric plane of the test mannequin were treated as replicate samples, yielding a total of 636 measurement locations to validate the CFD model.

Particle experiments were conducted using a Lechler ultrasonic spray system [US 1 nozzle (710.070.16.50 [EC] ) with a US 1 generator (071.091.01.11 [EC] ), Lechler, Inc., Metzengin, Germany], positioned horizontally upstream of the mannequin mouth at (x, y, z) of (–1.523, 0, 0.3685 m). The size distribution of Inland Oil (density = 0.86 g cc–1) aerosol from this ultrasonic nozzle was determined by optical sizing and counting in a settling chamber; the size-specific concentration inhaled into the mouth was also determined by optical sizing. The experimental aspiration fraction was calculated using the ratio: inhaled particle count per ml liquid aerosolized/generated particle count per ml liquid aerosolized, for a given particle diameter. Aspiration fractions were computed for particle sizes of 51.5–63.5 µm. These large particles were selected because the ultimate use of this work is to investigate the aspiration efficiency curve in low-velocity environments, and these particles are in the large scale of the inhalable spectrum.

Numerical
Gambit 2.1.6 (Fluent, Inc., Lebanon, NH) was used to generate the computational domain. Vertices from a generic file of a 3-D female head formed the basis for the computer representation of this mannequin. Simplifications to this file included the removal of the fine details associated with the ears, eyes and lips. The vertex coordinates were scaled to match the computer model to the dimensions of the experimental mannequin. Additional adjustments were made to modify the relative positions of the nose, lips and chin to better match the experimental mannequin. To reduce the size of the computational domain, lateral symmetry was assumed. Hence, a plane bisected the mannequin mid-sagittally, and only the right side of the mannequin head/neck was used for the numerical study. A total of 339 vertices were used to define the final head and neck geometry. The right half of a cylinder, bisected mid-sagittally, was dimensioned to the test mannequin to form the rest of the computational human form (Fig. 1).

The width of the computational domain was 0.762 m to match the half-width of the experimental wind tunnel. The height of the computational domain was 0.8212 m. The distance from the top of the domain to the top of the mannequin head was matched, but the torso was truncated at approximately hip level (0.25 m below the mouth opening) to save computational costs. The domain was 3.050 m long, with the mannequin mouth opening located 1.847 m into the domain, to match experimental particle inhalation conditions. The location of the mannequin was more than 10 head diameters upstream of the wind tunnel exit, sufficient to assume a stress-free boundary at the wind tunnel exit.

Because of the complex geometry associated with the anatomical shape in the wind tunnel, a tetrahedral meshing scheme with paved meshes was necessary. The spacing of nodes along the long edges of the wind tunnel domain was denser near the center of the wind tunnel, where the mannequin was located, and less dense near the wind tunnel entrance and exit to resolve velocity gradients more effectively. Three successively denser meshes were made by increasing the number of nodes on each edge by a factor of ~1.3, yielding a mesh-coarsening factor, r, of ~1.24. The three meshes contained ~159 900 (coarse), 308 300 (mid) and 591 300 (fine) nodes. Figure 2 illustrates the computational domain of the coarsest mesh.


Figure 2
View larger version (63K):
[in this window]
[in a new window]
 
Fig. 2. Illustration of coarsest mesh, N = 159 900 (a) domain, (b) close-up of face.

 
The equations used to model the airflow were steady state, incompressible, turbulent, Navier–Stokes equations:

Formula 1(1)

Formula 2(2)
The Boussinesq constitutive relationship was used to model the Reynolds stresses; the eddy viscosity was defined by:

Formula 3(3)
The standard two-equation k-epsilon model equations were used:

Formula 4(4)

Formula 5(5)
where k is the turbulence kinetic energy, {varepsilon} is the turbulent dissipation, {Phi} is the viscous dissipation function, and the standard constants are:

Formula 5
Equations 1, 2, 4, and 5 are the six equations used to solve the six degrees of freedom for the flow field: U (streamwise velocity), V (lateral velocity), W (vertical velocity), P (pressure), k (turbulence kinetic energy) and {varepsilon} (turbulent dissipation).

Table 2 provides a summary of the boundary conditions assigned to these simulations. On the domain inlet, a uniform velocity of 0.30 m s–1 with turbulence kinetic energy of 3.2912 x 10–5 m2 s–2 and dissipation of 3.58183 x 10–7 m2 s–3 were specified. These values were based on the experimental measurements of turbulence intensity (1.55%) and the assumption of 25 as Ru, the ratio of eddy viscosity to laminar viscosity. This Ru was slightly higher than the typical value of 10 used in wind tunnel studies to minimize the cell Reynolds numbers associated with these small turbulence intensities. The outlet was specified as no-stress, the lateral wall and top of the domain were specified as no-slip walls and the symmetry plane that bisected the mannequin laterally was specified to have zero normal velocity and zero tangential stress. Because the bottom of the computational domain was not at the same height as the floor of the wind tunnel, it was specified as having zero normal velocity and zero tangential stress to enhance the model's ability to produce the upward velocity as air approached the bluff body.


View this table:
[in this window]
[in a new window]
 
Table 2. Boundary conditions for numerical simulation

 
The mouth inlet was specified with normal velocity to achieve 10 lpm inhalation across the mouth inlet area and zero lateral and vertical velocity in that plane. For all other surfaces that defined the mannequin, the no-slip wall condition was assigned. Normal derivatives of turbulence kinetic energy (k) and turbulent dissipation ({varepsilon}) were set to zero at the mouth inlet. Because mesh refinement made slight adjustments to the mouth inlet size, minor modifications were made to the assigned normal velocity at the mouth to ensure identical flow rates for each mesh.

The simulations were conducted using the segregated solver with pressure projection and an element-Reynolds number relaxation scheme. Streamline-upwinding was used to stabilize the convective terms, yielding a nearly second-order accurate method with no crosswind artificial dissipation. The iterative solvers for the linear systems were conjugate gradient squared and conjugate residual methods. The clip was set to 1012 to account for the range of turbulence kinetic energy and dissipation in this system. The steady-state simulations required a total of ~120 h of computational time using a 3 GHz, 2 GB RAM computer (Windows 2000). Solutions to coarser mesh simulations were used as initial conditions for finer mesh simulations.

The solutions to the fluid flow model were then used to conduct particle aspiration studies using both laminar and turbulent particle trajectory methods. Laminar particle trajectory calculations relied on the mean velocity field estimates from the fluid model, whereas turbulent particle trajectories included random component to account for turbulence in the air. Particles were released from a position matching the location of the experimental nozzle relative to the mouth entrance (–1.523 m upstream and 0.3685 m above the mouth). Although particles from the experimental nozzle were emitted from a circular annulus, the corresponding simulations released particles from a horizontal line extending from the nozzle center (Y = 0 m) to 0.004646 m, the distance from the experimental nozzle's center to the edge of the annulus. No attempt was made to include the ultrasonic nozzle itself in the computational domain; thus, the turbulent wake behind the nozzle was ignored. Because experiments were conducted with 2 lpm air through the nozzle, an initial release velocity of 1.802 m s–1 was assigned to the particles to match the experimental air velocity at the outlet. Time steps of 0.00065 s were used in particle transport calculations. Durations of 6.5 s allowed particles to travel from the nozzle to the mouth or to bypass the mannequin head, sufficient to determine whether a particle was inhaled or not. Sensitivity to both the release location (+5 mm vertically, –5 mm streamwise) and to initial velocity (±11%) was examined to determine if fairly large experimental errors would result in significantly different numerical predictions of particle aspiration. For turbulent particle simulations at a given particle size, counts of aspirated particles and total particles released were accumulated until the addition of more particles have insignificant (<0.5%) changes to the aspiration fraction. For most particle sizes, 2000 particle releases were adequate.


    EXPERIMENTAL FLUID FLOW RESULTS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
The experimental data (U, V, W and K) were summarized by a mean and standard deviation based on replicate measurements for each of the 636 measurement locations. For each location, 2–6 replicate measures were available (mean = 2.2); these data include actual replicate measurements as well as the reflected velocity measurements from the left-hand (–Y) to the right-hand (+Y) side of the mannequin. Although dependence between the right and left side measurements exists, the impact is expected to be minimal and was considered a useful trade-off to data acquisition. Despite data reflection, the low sample sizes at each measurement location required correction for the computed standard deviations (Sokal and Rohlf, 1981Go).

Table 3 presents the minimum and maximum values for each degree of freedom, along with the fraction of measurement locations that had coefficients of variation (cv) <0.1, 0.25 and 1.0. The streamwise velocity (U) and the vertical velocity (W) were more precise than lateral velocity (V) and kinetic energy (K). The lateral velocities were often near zero, which contributed to the larger cv. In addition, when centered on the nose tip, the mannequin head was 2 mm wider on the right (+Y) side compared with the left (–Y) side. Because the mouth opening was centered on the nose tip, this anthropometric discrepancy may account for greater variability in the replicate measures of lateral velocity.


View this table:
[in this window]
[in a new window]
 
Table 3. Ranges and fraction of measurements below the indicated coefficient of variation (cv)

 

    COMPUTATIONAL VERIFICATION
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
Verification studies were performed on three sequentially denser meshes. For each mesh, the non-linear iterative convergence limit was established. Next, mesh convergence was assessed to determine the independence of the solutions with respect to the fine mesh. Locations were investigated to assess monotonic convergence by examining local orders of accuracy; GCIs were calculated to estimate the numerical uncertainty for the validation studies.

Non-linear convergence
Using the traditional L2 error norm, defined in equation (6), the Fidap software monitored the non-linear convergence for all nodes in the computational domain.

Formula 6(6)
where the subscript denotes the solution iterationnumber.

Convergence for a simulation was achieved and a solution to the fluid flow domain was provided when the global solution error (GSE) norm for each of the 6 degrees of freedom over all computational nodes reached a specified level (default = 10–3). To evaluate the adequacy of this global solution tolerance, the GSE was reduced by one order of magnitude, and the fluid flow field was recomputed to achieve the lower GSE tolerance. The resulting flow data were then assessed to determine whether the reduction of this tolerance affected the fluid flow field in the area of interest. The comparison between the fluid flow field at two tolerance levels was made by calculating a local L2 error norm. For this study, two local areas were investigated. The ‘velocity measurement’ area matched the 636 locations where experimental measures were made. The ‘nozzle region’ area was defined by 2533 upstream locations to incorporate a region where particle transport would occur (X = 0 to –1.6 m, Y = 0 to 0.5 m, Z = –0.3 to 0.12 m).

The global solution tolerance was reduced successively until insignificant changes occurred in the local relative L2 error norms for these regions. Table 4 provides the relative L2 error norms for each degree of freedom between solutions at GSEs of 10–3 and 10–4 and GSEs of 10–4 and 10–5 for these two regions. For a given mesh density, the local relative L2 error norm decreased with reductions in GSEs, although minor oscillations for the finest mesh were identified. The local L2 error norms for the velocity measurement region were <5% for each degree of freedom, which was small relative to any mesh errors when the global solution tolerance was set to 10–5. The error norms for the velocity measurement locations were less than for the nozzle locations, where the mesh was less refined.


View this table:
[in this window]
[in a new window]
 
Table 4. Relative L2 error norms for non-linear convergence tolerance for (a) velocity measurement region and (b) nozzle region upstream of inhaling mannequin

 
Mesh convergence
After an iterative convergence tolerance limit was established at 10–5, a local mesh convergence ratio, R2, was calculated for each degree of freedom based on the absolute L2 error norm. This error norm was computed by comparing the differences in the three meshes, using the equation from Stern et al. (2001)Go:

Formula 7(7)
where Formula 7

The term {varepsilon}j,k is the difference between the coarser (j) and finer (k) mesh level values for a given degree of freedom.

Table 5 provides the mesh convergence ratios for the series of three meshes at the same two study locations within the calculation domain. For the velocity measurement region, near the head, the regional mesh convergence ratios for all degrees of freedom were between 0 and 1, indicating monotonic convergence for solutions at both 10–4 and 10–5 global tolerance limits. For the upstream nozzle locations, the mesh convergence is monotonic for all terms except pressure and dissipation, which indicates that the upstream meshes are too coarse to solve the equations accurately.


View this table:
[in this window]
[in a new window]
 
Table 5. Three-mesh convergence ratios, R2, at sequential global solution errors (GSEs), for velocity measurement and nozzle locations

 
Even though mesh convergence ratios for the local data indicate convergence near the head, locally divergent or oscillatory behavior may also occur at specific locations (Stern et al., 2001Go). Thus, local mesh convergence ratios were calculated at each of the locations previously grouped. Of the 3816 possible local convergence ratios from the velocity measurement location, 17% exhibited monotonic convergence, while 73% remained oscillatory and 10% were divergent. The minimum local |R1| was 8.5 x 10–5, indicating that no node was considered completely converged using the convergence criteria of |R1| < 10–16 (Stern et al., 2001Go). Of the points that exhibited monotonic convergence, 10% had orders of accuracy (p) within the theoretical range of 1–2, where:

Formula 8(8)
For the nozzle location, 68% indicated monotonic convergence.

The subjective, statistical approach presented by Flynn and Eisner (2004)Go was also used to investigate convergence by fitting regression equations to the data from sequential mesh solutions. The relationship between predicted streamwise velocity in the finest and middle mesh is shown in Fig. 3, while Table 6 provides regression data for all degrees of freedom. The reported significant digits on Table 6 are provided to illustrate the tight confidence intervals on the regressed slopes and intercepts. Figure 3 illustrates that the disagreement between the predicted velocities occurred mainly in the higher velocity measures, which occurred near the mouth inlet. With the sequential mesh refinement, the slopes, intercepts and coefficient of determination (R2) all approach the ideal values of 1, 0 and 1, respectively, with very tight confidence intervals. The greatest residual error occurred near the mannequin head with the dissipation term. These results indicate that little changes in velocity and turbulence kinetic energy were evident between the middle and finest meshes, which demonstrates good ‘mesh independence’.


Figure 3
View larger version (14K):
[in this window]
[in a new window]
 
Fig. 3. Comparison of numerically predicted streamwise velocity between middle and fine mesh at velocity measurement locations (GSE = 10–5).

 

View this table:
[in this window]
[in a new window]
 
Table 6. Verification regression statistics for (a) velocity measurement locations and (b) nozzle locations

 

    VALIDATION
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
Velocity validation
‘Validation’ of a computational model is achieved when the numerical predictions agree with experimental measurements. Because both numerical and experimental uncertainties exist, validation is achieved when the uncertainty intervals of both predicted and measured values overlap for all degrees of freedom. For a given position, the experimental uncertainty is defined by the standard deviation of the measured velocity, whereas the numerical uncertainty is defined by the computed GCI for the finest mesh (Roache, 1998Go). These GCIs are available only for locations where monotonic convergence is observed.

Table 7 provides information on the fraction of monotonically convergent locations where the simulated variable was within 1, 2 or 3 SDs of the measured velocity component. The lateral velocity (V) had the strongest agreement: 82% of the simulated values were within 3 SDs of the experimental data. This finding is attributable to the larger experimental variability discussed earlier, which contributed to more locations having overlapping confidence intervals. The streamwise velocity (U) near the mannequin face exhibited reasonable agreement with the experimental velocity measurements: of the monotonically converging locations, 54% of predicted U were within 3 SDs of the experimental data. Of the monotonically converging locations, only 12% of predicted W were within 3 SDs of the experimental data.


View this table:
[in this window]
[in a new window]
 
Table 7. Fraction of numerical predictions within the indicated experimental uncertainty ({sigma}m) for monotonically convergent locations, regardless of order of accuracy

 
The agreements between the measured and simulated velocity data are illustrated graphically in Figs 4Go6. In each plot, the dashed diagonal line indicates perfect agreement between measured and simulated velocities while the solid line indicates the least-squares fit between the numerical and experimental data, regardless of monotonicity (i.e. all 636 paired data). The error bars indicate ±1 SD for the experimental velocity measurements and ±1 GCI for the numerically simulated velocities. The simulated streamwise velocities over-estimated the experimental U at values >0.13 m s–1; the simulated lateral velocities over-estimated the experimental V > 0 m s–1; and the simulated vertical velocities under-estimated the experimental W for all but the very large vertical velocity measurements, i.e. >0.11 m s–1. Even though the simulated velocities tended to over-estimate or under-estimate experimental results, the trends were well captured with coefficients of determination (R2) of 0.91, 0.90 and 0.84 for U, V and W, respectively. When all velocity component data are combined and analyzed together, as shown in Fig. 7, it is apparent that the central portion of the curve, where the bulk of the data exists, illustrates agreement between measured and predicted velocities. Table 8 provides regression statistics for the individual and combined velocity data, indicating that the slope is 14% greater than perfect agreement and the intercept is –0.011 for the combined data. Again, the reported significant digits on Table 8 are provided to illustrate the tight confidence intervals on the regressed slopes and intercepts. The coefficient of determination is improved to 97% for the combined data, indicating a strong predictive relationship for the model. Although the model was not able to predict all velocities within experimental accuracy, the model predicted the trends in the velocity field well, thereby indicating the usefulness of the model (Lasher, 2001Go).


Figure 4
View larger version (19K):
[in this window]
[in a new window]
 
Fig. 4. Comparison of numerically predicted and experimentally measured streamwise velocity (U).

 

Figure 5
View larger version (19K):
[in this window]
[in a new window]
 
Fig. 5. Comparison of numerically predicted and experimentally measured lateral velocity (V).

 

Figure 6
View larger version (17K):
[in this window]
[in a new window]
 
Fig. 6. Comparison of numerically predicted and experimentally measured vertical velocity (W).

 

Figure 7
View larger version (14K):
[in this window]
[in a new window]
 
Fig. 7. Velocity comparison using all velocity components.

 

View this table:
[in this window]
[in a new window]
 
Table 8. Regression statistics for measured versus simulated velocity, for the equation: Numerical = slope (measured) + intercept

 

Figure 8
View larger version (24K):
[in this window]
[in a new window]
 
Fig. 8. Comparison of numerically predicted and experimentally measured kinetic energy (KE).

 
Possible explanations for disagreement between the modeled and measured velocity flow include: (i) inappropriate truncation at the bottom of the calculation domain; (ii) differences between the computational face and the experimental mannequin and (iii) errors positioning the LDA at mouth center, resulting in a systematic bias in the experimental measurements. The first error is probably a noteworthy contributor to the poor agreement among the vertical velocity data. The truncation of the mannequin at hip level and the subsequent assignment of zero vertical velocity at this bottom plane may have prevented the upward flow near the bluff body that occurred experimentally. If the differences in facial detail represent a significant contribution to velocity field differences, then a larger problem exists: this finding implies that the large variability in facial features of real people has a greater impact on the velocity field than other factors.

The third possible explanation, an LDA positioning error bias, was examined using the numerical model. Because the traverse system that moved the laser and optics through the experimental wind tunnel was capable of moving in 1 mm increments, the positioning sensitivity was investigated at half this distance. Therefore, the simulated flow field was examined at locations ±0.5 mm laterally, vertically and streamwise from each of the velocity measurement positions (e.g. x' = x + 0.0005 m, y' = y 0.0005 m, z' = z). For each of the 636 validation locations, 27 nearby locations were investigated numerically. No single velocity field from these 27 sets of data performed any better than the original velocity measurement set. For each of the 636 velocity measurement locations, a position-based standard deviation ({sigma}p) was calculated using the predicted velocity components for the 27 adjusted positions. The standard deviations associated with the 27 adjusted-position velocity averages were less than the numerical uncertainty defined by the GCI. Paired t-tests revealed these differences were significant for the vertical velocity field (P = 0.0008). These uncertainties were incorporated into the uncertainty analysis by adding 1{sigma}p and 3{sigma}p to the GCI uncertainty term. Results are incorporated into Table 7. The inclusion of positional uncertainty into the numerical prediction uncertainty provided only moderate improvements in the fraction of validated locations for the velocity field. Thus, positional bias is not sufficient to address the differences in predicted versus measured velocity fields.

Kinetic energy validation
Within the velocity measurement locations, the simulated turbulence kinetic energy (KE) was orders of magnitude larger than the experimental values (mean experimental = 4.2 x 10–6 m2 s–2, mean numerical = 2.6 x 10–3 m2 s–2), as illustrated in Fig. 8. Only four of the 158 monotonically converging locations were successfully validated: these four locations had numerical uncertainties that were greater than the predicted turbulence KE, so that their confidence intervals included zero. These are illustrated in Fig. 8 for data with large horizontal error bars. These large numerical uncertainties were sufficient to allow for overlapping of measured and predicted uncertainty intervals for KE.

The fraction of validated locations increased from 3 to 80% when including the additional 0.5 mm positioning uncertainty to turbulence kinetic energy estimates. Using positional uncertainty ({sigma}p) and numerical uncertainty (GCI), 126 of the 158 monotonically converging locations had uncertainties larger than the mean predicted KE for the given location. This again resulted in the inclusion of zero in the uncertainty intervals, sufficient to include the measured kinetic energy within the confidence intervals of numeric predictions, despite the orders of magnitude differences in predicted versus experimental kinetic energy. More distressing is the lack of correlation between the predicted and measured values of the turbulence kinetic energy (R2 = 0.0009). Simulations at higher turbulence intensity (3%, with Ru = 10) were also conducted for these mesh series, and similar poor correlations were found (R2 = 0.0006).

The discrepancy between measured and simulated turbulence kinetic energy is probably attributable to the low particle seeding in experimental velocity measurements (Fig. 8). A mere 10 Doppler bursts were recorded per second, resulting in a filter cut-off frequency of 1.5 Hz, using the equation

Formula 9(9)
For this low particle-seeding rate, the true kinetic energy has been under-reported by LDA (Noback et al., 1998).

As a result of the experimental uncertainties, the turbulence kinetic energy from the CFD simulations remains unvalidated. Concerns of the thickness of the boundary layer associated with the high Reynolds number model are unresolved: in this low Reynolds number regime, the facial features of the human form may have been contained within the relatively thick viscous boundary layer present in low Reynolds number regimes, making resolution of kinetic energy near the walls difficult. Initial simulations using the low Reynolds number k-omega turbulence equations resulted in unachievable convergence at these same mesh densities. Finer near-wall meshes are required, and sensitivities to freestream values of turbulent frequency ({omega}) have not been resolved at this time.

Particle aspiration validation
Simulated particle releases occurred at locations representing experimental conditions: 1.523 m upstream and 0.3685 m above the mouth center, along a line from the plane of symmetry (Y = 0 m) to 0.004636 m laterally. Again, this work did not treat the aerosol as a uniform concentration within the wind tunnel, but rather used an aerosol source that was positioned to enhance the aspiration of larger particles in this small-scale system. Table 9 provides the aspiration fraction for experimental measurements as well as both laminar and turbulent particle aspiration simulations. For laminar particle simulations, no particles smaller than 55.0 µm and no particle larger than 58.0 µm released along this line were aspirated. Particles smaller than 55.0 µm hit the mannequin above the mouth and particles larger than 58.0 µm settled below the region of the mouth with the release location matching the experimental particle generation. This range was smaller than both the median aspirated particle size (60.2 µm) and the particle size with the maximum experimental aspiration fraction (66.5 µm). As indicated in Table 9, all of the 57.5 µm particles were aspirated by the simulated human.


View this table:
[in this window]
[in a new window]
 
Table 9. Particle aspiration fraction comparison

 
To compare numerical aspiration fractions with experimental results, turbulent particle motion was examined for 51.5, 54.5, 57.5, 60.5 and 63.5 µm oil particles (density = 0.86 g cc–1). These sizes matched binned data from optical particle count concentrations from experimental studies. Particles outside this range had negligible simulated aspiration fractions. The maximum simulated aspiration was 11.8% for 57.5 µm, within the range of the particle sizes inhaled in laminar simulations. The peak aspiration fraction for the experimental study occurred for 66.5 µm particles. Sensitivity tests indicated that neither moving the source upward 5 mm nor changing the initial particle velocity magnitude would shift the simulated aspiration fraction curve toward larger particles. The particle size differences between the measured and simulated aspiration curves are consistent with the under-prediction of vertical velocity identified earlier. The differences in the aspirated particle sizes indicate that differences between the modeled and measured velocity fields limit the model's usefulness to study aspiration efficiency from a point source.

Table 9 also illustrates that the simulated aspiration fraction was nearly two orders of magnitude greater than those observed experimentally. Poorly modeled turbulence kinetic energy and the oversimplification of the nozzle in the particle simulation study contributed to this finding. Because the fluid flow model over-estimated kinetic energy near the head of the mannequin, turbulent particle motion is expected to be more random and should result in fewer inhaled particles. However, this outcome is the opposite of the trend observed, and the KE discrepancies in the model are not the cause of these over-estimations in aspiration.

Therefore, the major contribution to the poor agreement between the magnitudes of aspiration fractions is probably the oversimplification of the particle generation mechanism. The source of air exiting the nozzle was represented only with the assignment of initial velocity to the particle, under-representing particle momentum. In addition, the lack of a realistic physical representation of the nozzle in the fluid flow field ignored the effect of the wake turbulence downstream of the nozzle. This wake effect increased particle dispersion experimentally, resulting in decreased particle aspiration at the downstream inhaling mannequin. Re-solving the fluid flow field would be necessary to fully capture this event. However, while better representation of the nozzle would probably yield better agreement between the magnitudes of the aspiration efficiency, it is not anticipated to produce better agreement in the shape of the aspiration fraction curves. The maximum aspirated particle size would still be smaller for the model than what was found experimentally.


    CONCLUSIONS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
The simplifications in the CFD model may have prevented better agreement between the experimental and simulated data. First, and perhaps most important, the truncation of the mannequin at hip height may have led to a consistent under-estimation of vertical velocity. The bottom of the computational domain presented here restricted the flow to zero-vertical velocity, allowing for streamwise and lateral velocity, unlike the no-slip condition of an actual wind tunnel floor. However, the location of this plane may have been insufficient to represent the upward vertical velocity field found in the experimental studies with a non-truncated mannequin. As a result of the reduced upward velocity, larger particles would settle faster toward the floor, resulting in a shift to slightly smaller particles being aspirated from the same point source. This is consistent with the 9 µm-shift in aspirated particle sizes identified in this study.

Previous work in particle inhalability proved the importance of the torso to induce upward particle flow in front of an inhaling mannequin (Vincent and Mark, 1982Go). Results from the present study indicate that the truncation of the mannequin at hip height may not adequately represent the velocity field associated with the full-height mannequin. As such, particle motion through this field may be affected by this reduction in upward velocity, even for the relatively large particles and low freestream velocity. Fortunately, the traditional method to determine numerical and experimental aspiration efficiency relies on the assumption of uniformly distributed particles within the wind tunnel. In this case, the sensitivity of the source position is avoided while still allowing examination of particle inhalability.

The second over-simplification was the particle release simulations, which led to over-estimation of the fraction of particles aspirated by the human form downstream of the ultrasonic nozzle. Additional modeling efforts to represent the aerosol source from this nozzle are necessary to validate the particle transport portion of this research. Fortunately, the particle sizes that were aspirated in the simulations were in the size range of those aspirated experimentally. The 9-µm shift to smaller aspirated particle sizes was consistent with the decreased vertical velocity field identified earlier.

The final simplification of the complex reality associated with particle inhalation studies was the selection of the computational mannequin shape. The computational human shape was a generic form, modified to match several key dimensions of the experimental mannequin, including lip location, mouth opening, head height and width, and location of shoulders. However, the specific facial details were not addressed, and, in light of the results discussed above, the effects of using the generic face shape could not be specifically addressed.

The greatest concern for this model is the lack of agreement between turbulence kinetic energies, but this is likely attributable to the measurement errors associated with low particle seeding. Validation is incomplete for turbulence kinetic energy simulations at this time.

Finally, this work provides recommendations for future CFD investigations that model airflow around an inhaling mannequin. First, the difficulties of obtaining appropriate mesh densities near the mouth and in the freestream are unavoidable. The average element dimension was 0.4 mm at the mouth and was 5.81 mm at the domain inlet. Verification studies found that discretization near the mouth and head was sufficient, but because upstream verification was incomplete, future studies should examine smaller discretization for the freestream. This would require either reducing the computational domain's size or using software that can access more than the 2GB RAM. In addition, setting the proper clip in these simulations was the key to obtaining converged meshes. The default of 108 was insufficient for this flow regime, and a value of 1012 was needed to achieve convergence to 10–5 GSE; this information should save valuable time in future computational studies with a complex human form.

Future CFD studies to investigate particle aspiration from a point source should include the full mannequin height and additional discretization for sources more than 10 head diameters upstream. However, the model developed here may still be useful for investigations that target information on particle concentration where the exact location of particle source is not significant. Particularly, the aspiration efficiency studies could use this model as they rely on the assumption of uniform particle concentrations upstream of the inhaling human but do not rely on the specific particle source location. Therefore, a human-scale version of this model may still provide insights into inhalable particle studies.


    NOMENCLATURE
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 


View this table:
[in this window]
[in a new window]
 
 

    ACKNOWLEDGEMENTS
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 
This research was supported by grants from the National Institute of Environmental Health Sciences (P30ES10126) and the National Institute for Occupational Safety and Health, Centers for Disease Control (#1 R01 OH07363). Its contents are solely the responsibility of the authors and do not necessarily represent the official views of NIEHS and NIOSH.

Received May 3, 2005; in final form July 27, 2005


    REFERENCES
 TOP
 ABSTRACT
 INTRODUCTION
 METHODS
 EXPERIMENTAL FLUID FLOW RESULTS
 COMPUTATIONAL VERIFICATION
 VALIDATION
 CONCLUSIONS
 NOMENCLATURE
 ACKNOWLEDGEMENTS
 REFERENCES
 

Anthony TR, Flynn MR, Eisner A. (2005) Evaluation of facial features on particle inhalation. Ann Occup Hyg; 49: 179–93.[Abstract/Free Full Text]

Baldwin PEJ, Maynard AD. (1998) A survey of wind speeds in indoor workplaces. Ann Occup Hyg; 42: 303–13.[Abstract/Free Full Text]

Brohus H. (1997) CFD simulation of personal exposure to contaminant sources in ventilated rooms. Proceedings of Ventilation '97, The 5th International Symposium on Ventilation for Contaminant Control, Global Developments in Industrial Ventilation, Ottawa, Canada, 14–17 September 1997. Vol. 1. pp. 215–26.

Brohus H, Nielsen PV. (1996) CFD models of persons evaluated by full-scale wind channel experiments. Proceedings of Roomvent '96, 5th International Conference on Air Distribution in Rooms, Yokohama, Japan, 17–19 July 1996. Vol. 2. pp. 137–44.

Chung IP, Dunn-Rankin D. (1992) Numerical simulation of two-dimensional blunt body sampling in viscous flow. J Aerosol Sci; 23: 217–32.

Dunnett SJ. (1997) A numerical study of the flow field in the vicinity of a bluff body with aspiration oriented to the flow. Atmos Environ; 31: 3745–52.

Dunnett SJ. (1999) An analytical investigation into the nature of the airflow near a spherical bluff body with suction. J Aerosol Sci; 30: 163–71.[Medline]

Dunnett SJ, Ingham DB. (1986) A mathematical theory to two-dimensional blunt body sampling. J Aerosol Sci; 17: 839–53.

Dunnett SJ, Ingham DB. (1988) The human head as a blunt aerosol sampler. J Aerosol Sci; 19: 365–80.

Dunnett SJ, Vincent JH. (2000) A mathematical study of aerosol sampling by an idealised blunt sampler oriented at an angle to the wind: the role of gravity. J Aerosol Sci; 31: 1187–203.

Erdal S, Esmen NA. (1995) Human head model as an aerosol sampler: calculation of aspiration efficiencies for coarse particles using an idealized human head model facing the wind. J Aerosol Sci; 26: 253–72.

Flynn MR, Eisner AD. (2004) Verification and validation studies of the time-averaged velocity field in the very near-wake of a finite elliptical cylinder. Fluid Dyn Res; 34: 273–88.[CrossRef]

Hayashi T, Yoshiaki I, Kato KS et al. (2002) CFD analysis of characteristics of contaminated indoor air ventilation and its application in the evaluation of the effects of contaminant inhalation by a human occupant. Building and Env; 37: 219–30.[CrossRef]

Hyun S, Kleinstreuer C. (2001) Numerical simulation of mixed convection heat and mass transfer in a human inhalation test chamber. Int J Heat Mass Transfer; 44: 2247–60.[CrossRef]

Ingham DB. (1981) The entrance of airborne particles into a blunt sampling head. J Aerosol Sci; 12: 541–9.

Ingham DB, Hildyard ML. (1991) The fluid-flow into a blunt aerosol sampler oriented at an angle to the oncoming flow. J Aerosol Sci; 22: 235–52.

Lasher WC. (2001) Computation of two-dimensional blocked flow normal to a flat plate. J Wind Engin; 89: 493–513.

Murakami S, Zeng J, Hayashi T. (1999) CFD analysis of wind environment around a human body. J Wind Engin; 83: 393–408.

Nobach H, Müller E, Tropea C. (1998) Efficient estimation of power spectral density from laser Doppler anemometry. Experiments in Fluids; 24: 499–509.[CrossRef]

Ogden TL, Birkett JL. (1977) The human head as a dust sampler. In Walton, W.H. editor. Inhaled particles IV: proceedings of an international symposium organized by the British Occupational Hygiene Society. Oxford: Pergamon Press. pp. 93–105. ISBN 0080205607

Roache PJ. (1998) Validation in Computational Science and Engineering. Albuquerque, NM: Hermosa Publishers.

Sokal RR, Rohlf FJ. (1981) Biometry. 2nd edn. San Francisco, CA: W.H. Freeman and Co.

Stern F, Wilson RV, Coleman HW et al. (2001) Comprehensive approach to verification and validation of CFD simulations—part 1: methodology and procedures. J Fluids Engin; 123: 793–802.[CrossRef]

Vincent JH, Mark D. (1982) Applications of blunt sampler theory to the definition and measurement of inhalable dust. Ann Occup Hyg; 26: 3–19.[Free Full Text]


Add to CiteULike CiteULike   Add to Connotea Connotea   Add to Del.icio.us Del.icio.us    What's this?



This Article
Right arrow Abstract Freely available
Right arrow FREE Full Text (PDF) Freely available
Right arrow All Versions of this Article:
50/2/157    most recent
mei048v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Add to My Personal Archive
Right arrow Download to citation manager
Right arrowRequest Permissions
Google Scholar
Right arrow Articles by ANTHONY, T. R.
Right arrow Articles by FLYNN, M. R.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by ANTHONY, T. R.
Right arrow Articles by FLYNN, M. R.
Social Bookmarking
 Add to CiteULike   Add to Connotea   Add to Del.icio.us  
What's this?