Characterization of Submicrometer Aerosol Deposition in Extrathoracic Airways during Nasal Exhalation

Abstract

Submicrometer and especially fine aerosols that enter the respiratory tract are largely exhaled. However, the deposition of these aerosols under expiratory conditions is not well characterized. In this study, expiratory deposition patterns of both ultrafine (<100 nm) and fine (100–1000 nm) respiratory aerosols were numerically modeled in a realistic nasal-laryngeal airway geometry. Particle sizes ranging from 1 through 1000 nm and exhalation flow rates from 4 through 45 L/min were considered. Under these conditions, turbulence only appeared significant in the laryngeal and pharyngeal regions, whereas the nasal passages were primarily in the laminar regime. Exhaled particles were simulated with both a continuous-phase drift flux velocity correction (DF-VC) model and a discrete Lagrangian tracking approach. For the deposition of ultrafine particles, both models provided a good match to existing experimental values, and simulation results corroborated an existing in vivo–based diffusion parameter (i.e., D ^0.5 Q ^−0.28). For fine particles, inertia-based deposition was found to have a greater dependence on the Reynolds number than on the Stokes number (i.e., St^0.1 _kRe^0.9), indicating that secondary flows may significantly influence aerosol deposition in the nasal-laryngeal geometry. A new correlation was proposed for deposition in the extrathoracic airways that is applicable for both ultrafine and fine aerosols over a broad range of nasal exhalation conditions. Results of this study indicate that physical realism of the airway model is crucial in determining particle behavior and fate and that the laryngeal and pharyngeal regions should be retained in future studies of expiratory deposition in the nasal region.

INTRODUCTION

Exposure to environmental aerosols of submicrometer size may represent a significant risk to human health. Submicrometer aerosols include particles in the ultrafine (<100 nm) and fine (100 nm to 1 μm) regimes. Recent studies indicate that aerosols in this size range are more biologically active and potentially more toxic than micrometer particles of the same material (Oberdorster and Utell 2002; Li et al. 2003; Kreyling et al. 2004; Kreyling et al. 2006). Sources of submicrometer aerosols include diesel exhaust (50 to 500 nm) (Kittelson 1998), cigarette smoke (140 to 500 nm) (Keith 1982; Bernstein 2004), and radioactive decay (1 to 200 nm) (ICRP 1994). Submicrometer bioaerosols include respiratory specific viruses such as Avian flu and SARS, which typically range from 20 to 200 nm (Mandell et al. 2004). These aerosols may deposit in the respiratory airways in discrete amounts resulting in local injury and the spread of infectious diseases.

Considering the extrathoracic nasal airways, which include the nasal passages, pharynx, and larynx, the deposition of submicrometer aerosols is associated with a number of detrimental health effects. The deposition of cigarette smoke particles has been quantitatively linked to the formation of respiratory tract tumors at specific sites (Martonen 1986). Yang et al. (1989) reported that respiratory tract cancers per unit surface area are approximately 3,000 times more likely in the extrathoracic airways including the larynx. Therefore, it is expected that cigarette smoke particle deposition per unit surface area is high in the extrathoracic region and contributes to cancer formation. A previous study by Cheng et al. (1993) has shown that ultrafine thoron progeny deposits at significant levels in the extrathoracic airways, which is associated with radiation induced tumor formation (ICRP 1994). Furthermore, a recent report by the National Research Council (NRC 2008), indicates that the in vivo response to inhaled bioaerosols is largely dependent on the site of local particle deposition within the respiratory tract. As a result of these potential health effects, the extrathoracic deposition of ultrafine and fine respiratory aerosols has been previously characterized during nasal inhalation (Cheng et al. 1988; Cheng et al. 1996b; Yu et al. 1998; Shi et al. 2006; Zamankhan et al. 2006; Liu et al. 2007; Xi and Longest 2008c). However, not all submicrometer aerosols deposit during inhalation. These remaining particles may then be deposited during exhalation throughout the respiratory tract, including within the extrathoracic region. Furthermore, deposition characteristics during exhalation are expected to be significantly different from inhalation conditions (Longest and Vinchurkar 2009).

A limited number of in vivo studies have considered the nasal deposition of submicrometer aerosols during exhalation in human volunteers. The study of Cheng et al. (1996a) considered the deposition of ultrafine aerosols ranging from 4–150 nm in 10 subjects and quantified variability in nasal geometries using MRI scans. Deposition during exhalation was shown to be a function of nasal cavity surface area, minimum cross-sectional area, and shape complexity. A correlation for nasal expiratory deposition was developed as a function of geometric parameters and particle diffusion. As a result, it was concluded that nasal deposition of ultrafine aerosols was highly variable among subjects and depended largely on the diffusional transport mechanism.

A common disadvantage of in vivo experimental methods for nasal deposition studies is difficulty in determining local deposition values. Imaging methods have become available that can resolve general regions of particle deposition within the nasal cavity (Lee and Berridge 2002). However, these methods have not been applied to the deposition of submicrometer aerosols in the nasal airways. Furthermore, no experimental in vivo data is available for the deposition of fine respiratory aerosols ranging from approximately 150 nm to 1 μm in the nasal cavity.

In vitro studies of submicrometer aerosol transport and deposition provide the advantage of avoiding human subject testing and can be used to determine deposition within general regions of the extrathoracic airways. Experimental studies that have evaluated expiratory deposition of ultrafine aerosols in replicas of the human nasal cavity include Cheng et al. (1995), Cheng et al. (1993), and Yamada et al. (1988). The nasal geometries used in these studies are typically derived from medical scan data (e.g., MRI) or casts of cadavers and extend from the nostrils to the upper trachea. Results of the available in vitro experiments are in general agreement with the deposition data from in vivo studies for ultrafine aerosols (Cheng et al. 1996a); however, less dependence on flow rate is often observed in the in vitro measurements (Cheng et al. 1996b). Cheng (2003) incorporated in vivo and in vitro nasal deposition data of ultrafine and micrometer particles to develop a correlation of expiratory nasal deposition for particle sizes in the diffusional and the impaction deposition regimes. In a series of studies, Kelly et al. (2004a; 2004b) evaluated the effects of surface roughness in nasal replica models that extended from the nostrils to the posterior pharynx and found a negligible influence on inspiratory deposition for particles in the diffusional regime. As with in vivo experiments, current in vitro models have also generally neglected the deposition of fine aerosols in the range of 200 nm through approximately 1 μm. Furthermore, extrathoracic aerosol deposition on a localized basis during exhalation has not previously been reported.

A number of numerical studies have considered the transport and deposition of fine and ultrafine particles in the nasal airways (Yu et al. 1998; Schroeter et al. 2001; Martonen et al. 2003; Shi et al. 2006; Zamankhan et al. 2006; Liu et al. 2007, Xi and Longest 2008c). Similar computational fluid dynamics (CFD) studies have also evaluated the transport and absorption of vapors in the nasal passages (Scherer et al. 1994; Hubal et al. 1996; Zhao et al. 2004). Comparisons of CFD results to experimental deposition data for submicrometer aerosols in the nasal cavity are often difficult due to differences in the geometric models and the complexity of the transport dynamics. Martonen et al. (2003) developed a numerical model of ultrafine aerosol transport and deposition in the nasal cavity based on diffusional theory in a highly simplified tubular geometry. The resulting correlation agreed reasonably well with experimental predictions of total deposition for nanoparticles ranging from 1–100 nm. Shi et al. (2006) reported local deposition patterns of 1 and 2 nm particles in the nasal cavity under laminar and transient breathing conditions. For submicrometer aerosols, Longest and Xi (2007a) showed that inertia may significantly affect the total deposition of particles as small as 70 nm and the local deposition of 40 nm particles in a bifurcating airway model. Moreover, the deposition rate in this diffusional-inertial regime is low (i.e., <1 ∼ 2%), rendering comparisons to experiments a more challenging computational problem. Zamankhan et al. (2006) simulated the inspiratory transport and deposition of ultrafine aerosols in the range of 1–100 nm and showed good agreement with the in vivo deposition data of Cheng et al. (1996a) for particles less than approximately 20 nm. However, particles greater than 20 nm did not match the experimental data to a high degree. Liu et al. (2007) simulated particles ranging from 0.354–16 μm using three Lagrangian particle tracking models. Results of this study showed that the standard eddy interaction model used for Lagrangian particle tracking in turbulent flows significantly over predicted the deposition of the smaller particles considered. Xi and Longest (2008c) implemented a continuous-phase drift flux velocity correction (DF-VC) model in a nasal cavity for inhaled submicrometer particles and obtained good agreement with existing experimental data (i.e., 5–150 nm). However, no available numerical study has considered the nasal or extrathoracic deposition of submicrometer aerosols in the inertial regime under expiratory conditions.

During exhalation, the influence of the upstream laryngeal and pharyngeal regions on transport and deposition in the nasal cavity may be significant. Longest and Vinchurkar (2007) showed the necessity of including upstream effects for validating CFD predictions of aerosol deposition in a multi-bifurcating airway geometry. Xi and Longest (2007b) evaluated the effects of physical realism on deposition patterns for micrometer aerosols in a mouth–throat model with varying degrees of geometric complexity. This study reported that geometric realism had a major effect on local inertia-based depositions and highlighted the importance of a realistic glottal aperture and angled trachea on deposition localization and particle profiles entering the lungs. Similar effects of geometric realism and inclusion of the larynx on downstream deposition have been reported by Brouns et al. (2007), Lin et al. (2007), Xi and Longest (2008a), and Xi et al. (2008).

The objective of this study is to better characterize the deposition of submicrometer aerosols in an anatomically accurate nasal-laryngeal airway model during nasal exhalation on a regional and highly localized basis. Specifically, existing correlations for nasal airway deposition in the ultrafine range are extended to evaluate the deposition of submicrometer aerosols where both impaction and diffusion are present (20 nm–1 μm). The nasal–laryngeal model used in this study is constructed based on medical images of healthy human subjects. To characterize aerosol transport and deposition in the extrathoracic airways, a recently developed drift flux velocity correction (DF-VC) model is considered, which accounts for both inertial and diffusional effects on submicrometer aerosol transport and deposition. Deposition results are compared with a standard chemical species (CS) model, which neglects particle inertia. Differences between these two models are used to isolate and evaluate the inertia-associated deposition rates, which are needed to develop an empirical correlation for fine respiratory aerosols. Considering that the DF-VC model has only been evaluated in a limited number of studies, deposition results are also calculated using an extensively tested Lagrangian tracking approach (Longest and Xi 2007b; Xi and Longest 2007b; Xi et al. 2008). Comparison of the continuous-field and Lagrangian models serves to corroborate the results of the DF-VC approach and to establish an expected range of submicrometer aerosol deposition values in the extrathoracic airways. For validation purposes, results of both the Lagrangian and DF-VC models are compared with available in vivo and in vitro deposition data under expiratory conditions.

METHODS

The focus of this study is to characterize the regional and local deposition of submicrometer aerosols in the extrathoracic airways during nasal exhalation using both discrete and continuous-phase particle tracking approaches. To facilitate this goal, an anatomically realistic nasal–laryngeal airway model was constructed based on medical images. Expected flow regimes ranging from laminar through turbulent (i.e., 4–45 L/min) were considered for both ultrafine (1–100 nm) and fine (100–1000 nm) monodisperse particles (Table 1). Details of the geometry construction, boundary conditions, and particle transport models are described below.

TABLE 1 Test conditions for the discrete and continuous-phase particle tracking simulations

Parameter	Range
dp (nm)	1–1000 nm
	Ultrafine: 1, 2, 3, 5, 10, 20, 40
	Fine: 100, 200, 400, 600, 800, 1000
Exhalation flow rate, Qexhale (L/min)	4, 7.5, 10, 15, 20, 30, 45
Cunningham correction factor, Cc	221–1.15
Particle diffusivity, D (cm^2/s)	5.44 × 10^−2–2.84 × 10^−7
Schmidt number Sc	2.92–5.61 × 10^5
Inlet particle Stokes number, St k	1.68 × 10^−8–1.0 × 10^−3
Inlet flow Reynolds number, Re	368–3,302

Construction of the Nasal Airway Model

The nasal cavity is composed of two narrow passages that are separated by the nasal septum. Each nasal passage features three curved fin-like airway protrusions known as the superior, middle and inferior meatus (Figure 1). The osseous tissues beneath each meatus are termed the superior, middle and inferior turbinates, respectively, and form the lateral wall of the main passage (Figure 1). Proceeding further into the airways, the two passages merge at the distal end of the nasal cavity and form the beginning of the nasopharynx. The floor of the nasal cavity is formed by the hard palate, which is also the roof of the mouth. To construct the airway model, MRI scan images of the nasal cavity region for a healthy non-smoking 53 year old male (weight 73 kg and height 173 cm) were used in this study. The MRI images were provided by the Hamner Institutes for Health Sciences (Research Triangle Park, NC) and consisted of 72 coronal cross-sections spaced 1.5 mm apart spanning the nostrils to the nasopharynx. This dataset was originally reported by Guilmette et al. (1989) and has since been used in a number of nasal particle deposition experiments (Guilmette et al. 1994; Cheng et al. 1995; Kelly et al. 2004a; 2004b) and simulations (Schroeter et al. 2006; Shi et al. 2006). The multi-slice images were segmented in MIMICS (Materialise, Ann Arbor MI) according to the tissue-air contrast and were converted into a set of coronal poly-lines that defined the nasal airway. Based on these contours, an internal nasal surface geometry was constructed in Gambit 2.3 (Ansys, Inc.). Following the same procedure, a pharyngeal–laryngeal geometry was constructed based on CT images of a healthy 34 year old male and was connected to the nasal cavity model (Figure 1). The final surface geometry was imported into ANSYS ICEM 10 (Ansys, Inc.) as an IGES file for meshing. Due to the high complexity of the model geometry, an unstructured tetrahedral mesh was generated with high-resolution prismatic cells in the near-wall region (Figure 1b).

FIG. 1 Computational nasal-laryngeal airway model: (a) medical image based surface model and (b) ICEM-generated computational mesh composed of approximately two million unstructured tetrahedral elements and a very fine near-wall pentahedral grid. The nasopharynx is divided into two passages by the uvula. The forward protruding angle formed by the epiglottis as well as the pharyngeal sinuses on both sides of the lower pharynx were retained in the current airway geometry.

Morphometric dimensions of the two nasal passages are shown in Figure 2 in terms of coronal cross-sectional area and hydraulic diameters (D_h = 4A_i /P_i ) as a function of distance from the nostrils. Here A_i is the coronal cross-sectional area, and P_i is the perimeter of the cross-section. The dip in cross-sectional area at a distance of approximately 30 mm from the nostrils corresponds to the nasal valve, as noted by an arrow in the figure. One interesting feature about the nasal passages is that while the available flow area (i.e., coronal cross-sectional area) increases downstream, the effective flow area (which is proportional to the hydraulic diameter) progressively decreases (Figure 2). This inverse relationship results from a consistent increase in the coronal perimeter, which indicates an increasing complexity of the nasal airway from the nostrils through the nasal passage. From Figure 2, dimensions of the right and left passages are similar but not identical. This disparity may reflect the fact that the nasal cavity is a dynamic structure, with dimensions changing due to environmental and breathing conditions.

FIG. 2 Dimensions of the left and right nasal passages as a function of the distance from the nose tip based on coronal cross-sectional area and hydraulic diameter, D_h . The nasal value is marked with an arrow for the left and right passages.

The uvula, a projection of tissue suspended from the soft palate, was preserved in the airway model. Movement of the uvula alters the airway structure that connects the nasopharynx and oropharynx. In this study, the uvula rests on the back of the throat and partially obstructs the inferior nasopharynx, resulting in two passages. Below the uvula, the pharynx is characterized by narrow and flat air channels, which converge into the wedge-shaped glottis in the laryngeal region. In order to approximate the in vivo airway morphology that was captured in the images, anatomical details such as the epiglottis and the two pharyngeal sinus regions on either side of the larynx were also preserved (Figure 1a).

During exhalation, airflow enters the upper trachea, travels through the larynx, pharynx, and nasal turbinates, and exits the nostrils. By retaining the laryngeal–pharyngeal region, more physiologically realistic inlet flow conditions are provided to the nasal cavity. As a result, the flow field and particle deposition characteristics in the nasal cavity considered in this study will better represent in vivo conditions compared with some previous studies that excluded the larynx and pharynx. Therefore, inclusion of the laryngeal–pharyngeal region allows for the results of this computational study to be directly compared with previous in vivo nasal deposition data. Considering deposition within areas of the nasal–laryngeal geometry, subregions considered include the larynx, pharynx, nasopharynx (NP), and nasal cavity (Figure 1a). The nasal cavity can be further divided into the turbinate region (TR), olfactory region (OR), nasal valve region (VR), and vestibule (Figure 1a). The relative effect of particle inertia on deposition is expected to vary significantly within these structures based on flow heterogeneity.

Boundary Conditions

Steady exhalation was assumed for all simulations with uniform velocity profiles and particle distributions at the tracheal inlet (Figure 1a). Confluent airstreams in the daughter bifurcations during exhalation produce secondary flows with two or more pairs of counter-rotating vortices. This complex pattern of secondary motion has an augmented dispersion effect that is similar to turbulence and results in relatively blunt cross-sectional velocity profiles in the trachea during exhalation (Fresconi and Prasad 2007; Longest and Vinchurkar 2009). Therefore, a uniform distribution of velocity at the laryngeal inlet appears to be a reasonable first order approximation to in vivo conditions. Furthermore, it is important to note that the in vitro experiments of Cheng et al. (1993; 1995) did not include the effects of upstream lung bifurcations. Initial particle velocities were assumed to match the local fluid velocity. Atmospheric pressure conditions were assumed at the nasal outlets. The airway surface was assumed smooth and rigid with no-slip (u_wall = 0) and perfect absorption conditions. In the body, the extrathoracic airway is covered with a thin layer of mucus, which captures particles at initial contact and clears them to the throat or nasal vestibule by mucocilliary movement within a time period of 10 to 15 min. Mass diffusion and metabolism of deposited particles may occur within the mucus layer and may change the zero-concentration conditions at the wall. However, due to the slow speed of the mucocilliary movement compared with the intranasal airflow and relatively low deposition rates, the no-slip and perfect absorption conditions are reasonable approximations in this study.

Discrete and Continuous-Phase Transport Equations

The flow conditions considered in this study are assumed to be isothermal and incompressible. The mean inlet Reynolds number at the trachea varies from 368 to 3,302. The maximum Reynolds number based on the hydraulic diameter of the glottal aperture is approximately 8,037. The onset of turbulence has been reported to occur at much lower Reynolds numbers in the complex geometries of the respiratory tract compared with circular ducts (Chan et al. 1980). Therefore, laminar, transitional, and fully turbulent conditions in the nasal–laryngeal model are expected. To approximate these multiple flow regimes, the low Reynolds number (LRN) k-ω model was selected based on its ability to accurately predict pressure drop, time-averaged velocity profiles, and shear stress for transitional and turbulent flows. Moreover, the LRN k-ω model was shown to provide an accurate solution for laminar flow as the turbulent viscosity approaches zero (Wilcox 1998).

The transport and deposition of submicrometer particles are simulated using both Lagrangian and Eulerian models. Both approaches consider finite inertial effects that may be significant for submicrometer particle deposition. The aerosols evaluated in this study had a tracheal Stokes number (St _k = ρ _p d ² _p C_cU/18μ _c D_h ) range of 1.68 × 10⁻⁸ to 1.0 × 10⁻³ and were assumed to be dilute and to not influence the continuous-phase, i.e., one-way coupled particle motion. Here ρ _p and d_p are particle density and diameter, respectively, C_c is the Cunningham correction factor, U is the nasal inlet velocity, and μ _c is the dynamic viscosity of the continuous phase. In our previous studies, the Lagrangian tracking model enhanced with user-defined routines was shown to provide close agreement with experimental deposition data in upper respiratory airways for both submicrometer (Longest and Xi 2007b) and micrometer particles (Xi and Longest 2007a). The discrete Lagrangian transport equations can be expressed

where v_i and u_i are the components of the particle and local fluid velocity, respectively, and τ _p (i.e., ρ _p d ² _p /18μ _c ) is the characteristic time required for a particle to respond to changes in the flow field. The drag factor f, which represents the ratio of the drag coefficient C_D to Stokes drag, is based on the expression of Morsi and Alexander (1972). Non-continuum slip effects on the drag of submicrometer aerosols are accounted for using the Cunningham correction (C_c ) factor (Hinds 1999)

where λ is the mean-free path of air, which is assumed to be 65 nm. The effect of Brownian motion on particle trajectories is included as a separate force per unit mass term at each time-step based on the approach of Li and Ahmadi (1992). The influence of non-uniform fluctuations in the near-wall region is considered by implementing an anisotropic turbulence model proposed by Matida et al. (2004), which describes wall-normal turbulent velocities as

In this equation, ξ is a random number generated from a Gaussian probability density function, k is the turbulent kinetic energy, and f_v is a damping function component normal to the wall for values of y ⁺ less than approximately 40. Here y ⁺ is the dimensionless wall distance defined as y ⁺ = u∗y/ν where y is the normal distance to the wall, u∗ is the near-wall friction velocity, and τ _w is the wall shear stress.

Considering the continuous-phase particle models, the transport equation for the chemical species (CS) approach can be written as

where c represents the mass fraction of submicrometer particles, υ _T is the turbulent kinematic viscosity, Sc _T is the turbulent Schmidt number, taken to be 0.9, and D is the molecular diffusion coefficient determined by the Stokes-Einstein equation (Hinds 1999). Perfect absorption at the wall, i.e., c_wall = 0, was assumed.

A primary feature of the drift flux approach is the inclusion of particle inertia effects through the particle velocity (v_j ) in the convection term (Manninen et al. 1996; Longest and Oldham 2008), i.e.,

The particle velocity is determined based on the slip velocity (νs _j ) between the particle and airflow as (Manninen et al. 1996; Longest and Oldham 2008)

For a continuous-field solution and small particles, as considered in this study ( 1 μm), the particle slip velocity can be determined in terms of the fluid pressure gradient as (Manninen et al. 1996; Longest and Oldham 2008)

where ρ _m is the mixture density.

The standard drift flux model calculates the particle slip velocity in Equation (7) based on conditions at control-volume-center locations. Values at control volume centers are expected to provide a reasonable estimate of slip velocity through a majority of the flow field. However, slip velocity conditions at wall surfaces are expected to be significantly different from control-volume-center values due to high near-wall particle deceleration. To improve the performance of the standard drift flux approach, Longest and Oldham (2008) proposed a velocity correction based on a sub-grid Lagrangian particle solution. Instead of fully resolving near-wall particle inertia on a continuous basis with excessive control volumes, the DF-VC model employs an analytic solution of particle velocity between the near-wall control-volume-center location and the wall surface. For low particle Reynolds numbers, the individual Lagrangian transport terms become linear and separable (Fan and Zhu 1998). An analytical solution is then applied for wall-normal values of particle velocity and particle position. Further details of the DF-VC model are provided in Longest and Oldham (2008) and Xi and Longest (2008b).

Based on the continuous-field approach for particle transport, the local mass deposition as a result of both diffusional and inertial effects is expressed as (Longest and Oldham 2008)

In the above expression, v_n represents the wall-normal particle velocity scalar, which is only considered for motion toward the wall, and A_l is the local surface area of the control volume face. The second term on the right-hand-side in the above equation accounts for particle inertial impaction, which is absent in the standard chemical species model.

Deposition Factors

The regional deposition fraction (Df) is defined as the ratio of particles depositing within a region to the number of particles entering the system. For the continuous-field models considered in this study, the deposition fraction (Df) is calculated as the ratio of the summed local mass deposition rate to the inlet mass flow rate

where the local wall mass flow rates [mdot]_w,l are available from Equation (8), and the summation is performed over N local surface areas. The total inlet mass flow rate of particles is calculated as

In the above expression, c_in is the constant inlet mass fraction of particles, u_in is the particle inlet velocity, and A_trachea is the tracheal inlet area.

Localized deposition can be presented in terms of a deposition enhancement factor (DEF), which quantifies local deposition with respect to the total or regional deposition rate. A deposition enhancement factor, similar to the enhancement factor suggested by Balashazy et al. (1999), for local region l can be defined as

where the summation is performed over the region of interest, i.e., the upper airway model. In the expression above, A_l and A_r represent local and regional areas, respectively. For the continuous-phase deposition predictions, the local area was assumed to be the size of wall-adjacent control volume faces.

Numerical Method and Convergence Sensitivity Analysis

To solve the governing mass and momentum conservation equations in each of the cases considered, the CFD package Fluent 6 was employed. User-supplied Fortran and C programs were implemented for the calculation of initial particle profiles, particle deposition factors, grid convergence, and deposition enhancement factors. A specific set of user-defined functions was developed for implementation of the Lagrangian and DF-VC models, which included modules for anisotropic turbulent effects, near-wall velocity treatment, fluid-particle slip velocity calculations, the mass fraction advection term, and mass deposition (Longest and Xi 2007b; Xi and Longest 2007a; Longest and Oldham 2008; Xi and Longest 2008b). All transport equations were discretized to be at least second order accurate in space. A segregated implicit solver was employed to evaluate the resulting linear system of equations. This solver uses the Gauss-Seidel method in conjunction with an algebraic multigrid approach for improving the calculation performance on tetrahedral meshes.

A grid sensitivity analysis was conducted by testing the effects of different mesh densities with approximately 620,800, 1,140,400, 1,975,600, and 3,212,000 control volumes while keeping the near-wall cell height constant at 0.05 mm. Increasing grid resolution from 1,975,600 to 3,212,000 control volumes resulted in total deposition changes less than 1% for both the DF-VC and Lagrangian models. As a result, the final grid for reporting flow field and deposition conditions consisted of approximately 1,975,600 cells with a thin five-layer pentahedral grid in the near-wall region and a first near-wall cell height of 0.05 mm.

For discrete Lagrangian tracking, the number of seeded particles required to produce count-independent depositions was considered. Particle count sensitivity testing was performed by incrementally releasing groups of 10,000 particles. The number of groups was increased until the deposition rate changed by <1%. Due to the low deposition rates, more particles were required for fine aerosols to generate count-independent results compared with ultrafine aerosols. The final number of particles tracked for 1–40 nm and 100–1000 nm aerosols were 150,000 and 600,000, respectively.

RESULTS

Flow Fields

Figure 3 shows the steady state velocity field in the nasal-laryngeal geometry for an exhalation flow rate of 20 L/min. Air exhaled from the lungs enters the larynx, travels through the pharynx and nasal passages, and exits the nostrils. In addition to the geometric curvature from the larynx through the nasopharynx, the airflow experiences two dramatic geometric constrictions before it enters the nasal cavity, which arise from the glottis and the nasopharynx flow division formed by the uvula. Recirculation zones are observed both downstream of the glottis and within the nasopharynx. As a result, the inlet flow field entering the nasal cavity is far from uniform, with peak velocities located in the dorsal region of the nasopharynx and low flow near the center (see Slice 3-3′). As a result, peak airflow in the nasal cavity occurs in the middle meatus close to the septum, while a second peak occurs along the nasal floor of the inferior meatus. Only a minimal fraction of the total expiratory nasal airflow passes through the superior meatus where the olfactory epithelium is located.

FIG. 3 Velocity fields in the nasal-laryngeal airway at an exhalation flow rate of 20 L/min in three-dimensional and cross-sectional (coronal) views. The secondary motions within the nasal valves and middle turbinate region are also presented.

To illustrate the effect of secondary velocity motion, two-dimensional stream traces are shown in selected portions of two coronal slices (Slices 1-1′ and 2-2′ of Figure 3). Due to the thin air channels, vortices are largely damped within the nasal passages. However, a vortex is observed near the middle of the valve region, as illustrated in Slice 1-1′. The magnitude of the secondary velocity in each slice is approximately 30% of the axial flow. This secondary velocity component functions to distribute the exhaled air into the fin-like meatus passages. Of particular interest is the observation that these secondary velocities also distribute flow toward the olfactory region. As air progresses in and out of the projecting airways, the resulting secondary streams move upward prior to the olfactory region (Slice 2-2′) and downward after the olfactory region (Slice 1-1′) during exhalation. It is expected that this secondary motion near the olfactory region is delicately balanced in order to convey a sufficient amount of particles or vapors for perception while remaining small enough to protect this extremely sensitive area that is directly connected to the brain.

The evaluation of turbulence in the nasal airway during exhalation is shown in Figure 4 based on iso-surfaces of the turbulent viscosity ratio, ζ, which is the ratio of total to laminar viscosity ζ= (ν + ν _T )/ν. In this expression, ν and ν _T are the laminar and turbulent kinematic viscosities, respectively. In general, turbulence is observed mainly within the laryngeal and pharyngeal regions, and at a lower intensity within the posterior turbinate region. Relatively weak turbulence is observed in the main nasal passage. Figure 4a shows an iso-surface of ζ = 2, which implies equivalent laminar and eddy viscosities. This contour level is observed in the whole region of the larynx and pharynx, but extends through only part of the middle meatus in the turbinate region. Figure 4b displays conditions for which the turbulent viscosity is one order of magnitude above the laminar viscosity (ζ = 11). In this scenario, the iso-surface is limited primarily to the laryngeal–pharyngeal region, suggesting that the major portion of the nasal cavity can be regarded as laminar.

FIG. 4 Iso-surfaces of the turbulent viscosity ratio, ζ= (ν+ν _T )/ν, in the nasal-laryngeal airway at an exhalation flow rate of 20 L/min with (a) ζ= 2, and (b) ζ= 11.

Particle Transport

Discrete particle dynamics in the nasal airway are visualized in Figure 5 as both flow element streamtraces and snapshots of particle locations at selected instants during a washout phase at an exhalation flow rate of 20 L/min. In Figure 5a, six points were seeded at the inlet (i.e., upper trachea) and their trajectories were tracked through the airway geometry. Vortical flows are apparent in regions around the epiglottis, oropharynx, and dorsal nasopharynx. Due to enhanced turbulent mixing within the larynx and pharynx, fluid particles initiated on one side of the tracheal inlet (e.g., right side of the trachea in Figure 5a) exit through both sides of the nasal passages. In Figure 5b, three thousand 400 nm particles were released at t = 0 s, and particle positions after different periods of time were recorded. Rapid aerosol transport is observed in the laryngeal region due to flow constriction from the triangular-shaped glottic aperture, as illustrated at t = 0.01 s and 0.02 s. The particle profiles start to bifurcate at approximately t = 0.025 s because of the uvula obstruction and merge in the upper nasopharynx. Due to substantial flow reversals and recirculation in the upper nasopharynx, a portion of the exhaled aerosol will be retained, which may be breathed into the lungs during subsequent inhalations or may deposit due to diffusion and sedimentation. The seemingly random particle locations after t = 0.03 s indicate strong influences from the local vortical flow structures and possible turbulent dispersions. It is noted that the local particle distribution is a function of not just the local flow conditions but of the entire flow history the particles have experienced.

FIG. 5 Particle dynamics illustrated by (a) streamtraces of six fluid elements initiated at the upper trachea, and (b) snapshots of particle locations inside the nasal-laryngeal airway at different instants and a constant exhalation rate of 20 L/min. Three thousand discrete particles of 400 nm diameter were released at t = 0 s and tracked in the particle history illustration.

Particle Deposition

Figure 6 provides numerically determined deposition fractions based on the DF-VC and Lagrangian models in comparison to in vitro deposition data obtained in the nasal-laryngeal geometry during exhalation. The test conditions used in the in vitro experiments shown in Figure 6 are described in Table 2 where A_s is the nasal surface area and is the average cross-sectional area of the nasal passage. It is noted that the computational geometry in this study and the nasal component of the ANOT2 (Adult-Nasal-Oral-Tracheal) model were both constructed from the same set of nasal images. The in vitro and computational deposition data is for the extrathoracic region during exhalation including the larynx through the nasal cavity. As shown in Figure 6, both the DF-VC and Lagrangian methods appear to provide reasonably good agreement with in vitro measurements for the flow rates considered (i.e., 10 and 20 L/min). Specifically, for an expiratory flow rate of 10 L/min, the simulation results agree with the ANOT1 and ANOT2 data (Cheng et al. 1995) to a high degree. Furthermore, the close match between the continuous-phase DF-VC and discrete Lagrangian tracking methods indicates that with appropriate near-wall treatments, both models provide effective approaches in capturing the diffusional and inertial deposition considered. The slight underestimation of the simulation results compared with measurements may be attributed to particle electrostatic effects in the replica experiments, which enhance the particle deposition fraction (Cohen et al. 1998). Wall surface roughness in the experimental models may also increase deposition (Oldham 2006; Shi et al. 2007).

FIG. 6 Comparison of particle deposition fractions as a function of particle diameter based on the DF-VC and Lagrangian modeling approaches as well as existing in vitro experiments at an exhalation flow rate of (a) 10 L/min and (b) 20 L/min.

TABLE 2 Comparison of in vitro and in vivo experiments with the nasal-laryngeal geometry

		Original data	Measurement			Particle
Reference	Name	source/subject	protocol	(cm^2)	As (cm^2)	(nm)
In vitro Cheng Y.S. et al. (1993)	C1	Replica cast from cadaver specimen	Trachea in/nose out	NA	NA	1.2, 1.7
In vitro Cheng K.H. et al. (1995)	ANOT1∗	Replica cast from cadaver specimen	Trachea in/nose out	5.56	NA	3.6–100
In vitro Cheng K.H. et al. (1995)	ANOT2∗	Nasal cavity: MRI of 53 yr old male	Trachea in/nose out	4.80	NA	3.6–150
		Larynx-pharynx: cadaver specimen
In vivo Cheng Y.S. et al. (1996b)	NA	4 normal males (ages 36–57 yr)	Mouth in/nose out with oral tube	2.81–3.16	184–238	5–100
In vivo Cheng K.H. et al. (1996a)	NA	10 normal males (ages 24–58 yr)	Mouth in/nose out with oral tube	NA	190–264	4–150
Computational Current study	NA	Nasal cavity: MRI of 53 yr old male	Trachea in/nose out	3.13	237.5	1–1000
		Larynx-pharynx: CT 36 yr old male

To further evaluate the simulation results predicted with the DF-VC and Lagrangian models, airway deposition fractions from Figure 6 are plotted as a function of different existing diffusion parameters that have been suggested based on in vitro and in vivo deposition studies (Figure 7). Exhalation flow rates considered in Figure 7 include 4, 7.5, 10, 15, 20, 30, and 45 L/min with particle sizes ranging from 1 to 1000 nm. The extrathoracic region of deposition considered includes the larynx through the nasal cavity. Figures 7a and 7b show the variation of the DF-VC tracking results as a function of two diffusion parameters. The first parameter (D ^0.5 Q ^−0.125) was theoretically derived by Cheng et al. (1988) based on the assumption of turbulent diffusion in pipe flow and was later adopted in a series of in vitro studies (Cheng et al. 1993; Cheng et al. 1995; Cheng et al. 1996b) correlating diffusional deposition data within nasal airway replicas. The second parameter (D ^0.5 Q ^−0.28) was later derived by Cheng et al. (2003) based on in vivo nasal deposition data, which exhibited a greater dependence on flow rate (i.e., exponent of −0.28) than the replica-based parameter (i.e., exponent of −0.125). It is noted from Figure 7a that the DF-VC tracking results do not fully collapse for different breathing activities and particle sizes, suggesting that the replica-based quantity (D ^0.5 Q ^−0.125) does not fully account for the effect of exhalation flow rate. In contrast, a more precise correlation was obtained from plotting the simulation results as a function of the in vivo-based parameter (D ^0.5 Q ^−0.28). Following the format of the in vivo empirical correlation suggested by Cheng et al. (2003), Df = 1 − exp(-aD ^0.5 Q ^−0.28), a coefficient of a= 11.3 was obtained for the numerical data (Figure 7b). The resulting R ² value was 0.901, indicating reasonably good agreement between the numerical data and algebraic expression. However, noticeable deviations are still observed in the range of D ^0.5 Q ^−0.28= 2 × 10⁻³−6 × 10⁻². To potentially improve agreement with the DF-VC model results, the Cunningham correction factor C_c , which accounts for particle diameter and mean free path effects, was incorporated separately into the diffusion parameter as C^b _cD ^0.5 Q ^−0.28. Again, an equation of Df = 1 − exp(-cC^b _cD ^0.5 Q ^−0.28) was used to fit the expiratory simulation data. The best-fit values of b and c were 0.36 and 1.62, respectively, resulting in R ²=0.983 (Figure 7b). The resulting diffusion regime correlation for the deposition fraction of ultrafine particles in the nasal-laryngeal model during exhalation is expressed

As shown in Figure 7b, inclusion of the Cunningham correction factor results in a significantly better approximation of the numerical data. It is noted that this correlation is only recommended for ultrafine particles (diffusional regime), and will be further extended to the range of fine particles (diffusional-inertial regime) by including the inertial effects in a later part of this study.

FIG. 7 Comparison of the deposition fraction predicted by the DF-VC model with diffusional parameters based on existing (a) in vitro (D ^0.5 Q ^−0.125) and (b) in vivo (D ^0.5 Q ^−0.28) studies. Inclusion of the Cunningham correction factor (C_c ) improves agreement between the correlation and numerical data (b). Results of the DF-VC model and Lagragian particle tracking are compared with existing in vitro and in vivo deposition data for particles less than 150 nm as a function of (c) the in vivo diffusion parameter and (d) the diffusion parameter suggested by Cheng et al. (1996b). In general, deposition results from both the DF-VC and Lagrangian models match the experimental data. Units: Q [L/min] and D [cm²/s].

Figures 7c–d show the numerically determined deposition fractions using the DF-VC and Lagrangian models in comparison with available in vitro and in vivo deposition data, which are summarized in Table 2. Again, simulation results for the continuous-phase and discrete tracking approaches match closely. Reasonably good agreement is observed between the simulations and in vitro measurements for the three replica casts considered with slight under predictions from the numerical models, as shown in Figure 7c. Specifically, the closely correlated association of deposition fraction with the parameter D ^0.5 Q ^−0.28 indicates a stronger dependence of nasal deposition on particle diffusivity (exponent of 0.5) than on flow rate (exponent −0.28) for the submicrometer particles considered. Figure 7d shows the simulation results plotted based on a geometric deposition parameter and compares favorably with an empirical correlation from the in vivo nasal deposition study of Cheng et al. (1996b). Values of these geometric features are summarized in Table 2.

Finite Particle Inertia and Deposition

Figure 8 shows deposition fractions for the DF-VC and Lagrangian models in the nasal-laryngeal geometry as a function of the particle Stokes number (St _k ). In contrast to close agreement between these two models for diffusional regime particles, the Lagrangian approach generally predicts higher deposition for larger (diffusional-inertial regime) particles. Simulation results using a standard chemical species (CS) model are also reported in Figure 8 (“×” symbol) and can be used to estimate the purely diffusional component of deposition for fine aerosols. As described earlier, the DF-VC model accounts for both particle inertia and diffusion, whereas the CS model only considers diffusion. As a result, the increased deposition of the DF-VC model over the CS approach can be entirely attributed to particle inertia effects. For ultrafine particles where inertial effects are negligible, the predicted depositions appears identical for the two continuous-phase models considered, as expected. Noticeable deviations of deposition predictions between the DF-VC and CS models begin around St _k = 1.0 × 10⁻⁵ and become progressively larger with increasing St _k . Accordingly, a critical value above which inertia significantly influences deposition in the nasal cavity during exhalation is identified as St _k = 1.0 × 10⁻⁵. This limit is equivalent to a 90 nm particle under resting conditions (15 L/min) and a 50 nm particle under moderate activity conditions (30 L/min). Compared with the inertial limit of St _k = 5.0 × 10⁻⁵ for tracheobronchial airways (Longest and Xi 2007a), this smaller value in the nasal airway indicates an earlier onset of particle inertial effects, which may be due to the increased complexity of this geometry. Obviously, neglecting the finite inertia (i.e., as with the CS model) underestimates the deposition rate for particles with St _k > 10⁻⁵. It is interesting to note that the DF-VC deposition rate increases abruptly as the flow changes from laminar (Q = 4, 7.5 L/min) to transitional and turbulent (Q = 10∼45 L/min), suggesting that enhanced secondary flows and turbulent dispersion are potentially significant for fine aerosol deposition.

FIG. 8 Comparison of numerically determined deposition fractions as a function of Stokes number, St _k , for the DF-VC, Lagrangian, and CS models. Deviations in deposition fractions between the models that account for inertia (DF-VC and Lagrangian) and the chemical species (CS) approach occurs around St_k = 1.0 × 10⁻⁵ and indicates the onset of particle impaction.

While various forms of empirical deposition correlations exist for very small (<100 nm) and large (1 μm and greater) particles, correlations applicable for particles between these limits have rarely been reported. In this study, a correlation is sought that extends the diffusional expression of Equation (12) by incorporating inertia-associated effects. A significant issue in accounting for inertial effects is determining the appropriate parameter for correlating the data. Figure 9 shows the variation of inertia-based deposition rates plotted against two different parameters, i.e., (d ² _p Q) and (d ^0.2 _p Q). The inertial deposition rates were determined as the difference between the DF-VC and CS model predictions for the nasal-laryngeal geometry. A large scatter was found in Figure 9a for the conventional impaction parameter (d ² _p Q∼ St _k ) that has been extensively used to correlate microparticle deposition. However, a consistent pattern was noticed among the scattered data points with the inertial-based deposition values forming multiple distinct lines (Figure 9a). These multiple regular curves suggest that certain flow-rate-related influences are not fully considered. To better account for the exhalation flow rate, a new parameter with the form of (d^m _pQⁿ ) was tested using least-squares curve fitting. The parameter that was found to provide the best fit was d ^0.2 _p Q, which is directly proportional to St ^0.1 _k Re ^0.9. The improved data-fit based on this parameter, as shown in Figure 9b, indicates a much stronger dependence of deposition on Reynolds number (Re exponent of 0.9) than on inertial impaction (St _k exponent of 0.1) and may suggest that mechanisms other than diffusion and impaction need to be considered for correlations of fine respiratory aerosols. This result is somewhat expected considering that fine particles are largely influenced by flow features (i.e., secondary flows and turbulent dispersion).

FIG. 9 Inertia-based deposition fraction as a function of (a) the conventional impaction parameter (d ² _p Q), versus (b) a new inertial parameter (d ^0.2 _p Q) for particles ranging from 1 to 1000 nm and exhalation flow rates ranging from 4 to 45 L/min.

Based on the analysis above, a nasal-laryngeal airway deposition correlation applicable for both ultrafine and fine particles during exhalation is proposed as

and is used to correlate deposition values for particles in the range of 1 through 1000 nm and flow rates ranging from 4 through 45 L/min. Deposition values within these limits were determined using the DF-VC model. As shown in Figure 10, the best-fit values of coefficients a and b were 1.62 and 7.1 × 10⁻⁵, respectively, resulting in R ²= 0.94. The deposition of submicrometer aerosols in the nasal–laryngeal geometry during exhalation can then be expressed as

where the variable dimensions are as follows: diffusion coefficient D (cm²/s), exhalation flow rate Q (L/min), and particle diameter d_p (μm). Alternatively, this correlation can be expressed in terms of nondimensional variables as

where Sc = ν/D, Re = UD _h /ν, and St _k = ρ _p d ² _p C _c U/18μD _h are the Schmidt, Reynolds, and Stokes numbers, respectively. The characteristic length adopted in the above equation is the hydraulic diameter of the nasal valve, D_h = 4.7 mm as shown in Figure 2, and characteristic velocity was at the tracheal inlet.

FIG. 10 Deposition fraction as a function of the new composite deposition parameter over a range of flow rate and particle conditions.

Local and Sub-Regional Deposition Results

To highlight the influence of finite particle inertia on deposition localization, a comparison of deposition enhancement factors (DEF) predicted by the DF-VC model is shown in Figure 11 for four particle sizes at an exhalation flow rate of 20 L/min. As discussed, the DEF quantifies particle localization with respect to total deposition. Differences of the DEF distribution among these four particle sizes are significant. In contrast to the more uniformly distributed DEF values for 40 and 100 nm particles, DEF values for 400 and 1000 nm particles are highly heterogeneous and localized. A careful comparison between Figures 11(a) and 11(b) also reveals elevated localizations for 100 nm particles compared with 40 nm particles in certain areas (e.g., within the circled region, DEF = 15.6 for 100 nm vs. 6.4 for 40 nm). Furthermore, the maximum DEF values significantly increase as the particle size increases. For 400 and 1000 nm aerosols, highly localized depositions, or hot spots, were observed to closely associate with high-speed-flow or stagnation regions such as the glottis, epiglottis, and uvula. Specifically, around the nasopharynx division, i.e., where the uvula rests, the DEF values for 400 and 1000 nm aerosols are between one and two orders of magnitude larger than those for 40 nm particles. These hot spots may have important implications in assessing the etiology of epithelial lesions in the upper airways resulting from exposure to submicrometer toxic particles.

FIG. 11 Deposition enhancement factors (DEF) based on the DF-VC model for different particle sizes at an exhalation flow rate of 20 L/min: (a) 40 nm, (b) 100 nm, (c) 400 nm, and (d) 1,000 nm.

Partitionings of nasal airway deposition between the nose, pharynx, and larynx are shown in Figure 12 for three exhalation flow rates (10, 20, 30 L/min) and four particle sizes (40, 100, 400, and 1000 nm) based on the DF-VC model. Deposition percentages are also summarized in Table 3. The deposition ratio of each region is not constant, but varies with breathing conditions and aerosol sizes. As particle size increases from 40 to 1000 nm, diffusion is progressively decreased whereas inertia-associated effects increase. The minimum deposition fraction for 400 nm particles indicates that both diffusion and inertia are at a minimum for this size. It is interesting to compare the flow-deposition slopes at each region for the four particle sizes considered. For 40 nm particles, deposition decreases as the flow rates increases, which is expected due to diffusion-dominate wall contact. However, the nasal cavity shows the highest sensitivity (steepest slope) among the three regions, likely owing to slower moving flows (i.e., prolonged residence times) and large surface area in the nasal turbinate region. For 1000 nm particles, the relative inertial effects in the larynx and nasophaynx (included with the nose) are more pronounced than in the pharynx.

FIG. 12 Regional partitioning of deposition fractions between the nose, pharynx, and larynx under three expiratory conditions (10, 20, and 30 L/min) for different particle sizes: (a) 40 nm, (b) 100 nm, (c) 400 nm, and (d) 1,000 nm.

TABLE 3 Deposition fraction (%) in the nasal cavity, pharynx, and larynx during exhalation

Region	2.1. Q exhale (L/min)	2.2. Deposition fraction (%)
		40 nm	100 nm	400 nm	1,000 nm
Nasal cavity	10	1.193	0.392	0.089	0.134
	20	0.810	0.305	0.086	0.169
	30	0.622	0.274	0.091	0.197
	10	0.554	0.294	0.081	0.130
Pharynx	20	0.235	0.171	0.058	0.141
	30	0.240	0.146	0.053	0.148
	10	0.487	0.176	0.056	0.084
Larynx	20	0.412	0.212	0.059	0.110
	30	0.404	0.217	0.076	0.158

DISCUSSION

In this study, expiratory deposition in a scan-based nasal-laryngeal model was evaluated for fine and ultrafine aerosols using standard Lagrangian aerosol tracking and a new continuous-field drift flux approach. For ultrafine particles, a better correlation of predicted deposition was obtained with a diffusion parameter derived from in vivo results compared with a replica-based parameter, suggesting that the numerical modeling in this study captured major aspects of in vivo conditions. Furthermore, inclusion of the Cunningham correction factor in the diffusion parameter provided an improved fit to the numerical deposition data. For fine aerosols in the range of 100–1000 nm, deposition was found to be more dependent on flow rate than aerosol size. The resulting impaction parameter had the form d ^0.2 _p Q, which is proportional to St^0.1 _k Re^0.9. This relationship between flow rate and inertial deposition occurred for both laminar and turbulent flows, indicating that increased deposition as a function of flow rate was not only an effect of turbulent dispersion. Based on improved parameters in both the diffusional and impaction regimes, a new correlation was proposed to capture submicrometer aerosol deposition in the extrathoracic region during nasal exhalation over a broad range of flow rates. Numerical results also indicated significant particle localization with maximum deposition enhancement values ranging from 6.4 to approximately 50 for submicrometer aerosols.

Deposition results for ultrafine particles were found to correlate well with a diffusion parameter established from in vivo experiments (i.e., D ^0.5 Q ^−0.28). In contrast, a previously suggested diffusion parameter based on in vitro experiments did not result in a good correlation of the numerical data. One potential explanation of the improved correlation with the in vivo diffusion parameter is the inclusion of a scan-based larynx and pharynx. Flow features generated in the laryngeal and pharyngeal regions propagate into the nasal cavity during exhalation and influence aerosol transport and deposition. Several previous in vitro experiments on nasal deposition have included the laryngeal and pharyngeal regions (Yamada et al. 1988; Cheng et al. 1993; Cheng et al. 1995). However, it is not clear if the cast models used in the studies of Cheng et al. (1993, 1995) create flow disturbances consistent with the scan-based geometry used in this study. For example, the laryngeal and pharyngeal dimensions of the cadaver-derived ANOT1 model are expected to be substantially larger than those of the current scan-based model. Considering that the results of the current study correlate well with the in vivo diffusion parameter, the narrow airways of the scan-based geometry may be more consistent with conditions in living subjects.

A second potential reason for the difference in diffusion parameters observed between the numerical model and previous in vitro results may be based on geometric differences in the nasal cavity. As observed in Table 2, the average cross-sectional area () of the scan-based model used in this study is within the range of in vivo values reported by Cheng et al. (1996a). In contrast, the two replica casts, ANOT1 and ANOT2, employed by Cheng et al. (1995) had average cross-sectional areas that were 56% and 40% greater than the scan-based geometry, respectively, and not within the in vivo range of values reported by Cheng et al. (1996a). A similar difference was observed by Guilmette et al. (1989), who reported an approximate 40% increase in cross-sectional area for a cast-based model compared with MRI images from a living subject. The larger cast geometries may be a result of post-mortem shrinking of mucus tissue or expansion of the upper airway soft tissue due to the casting process. Based on these observations, it appears that the scan-based model may be geometrically more similar to in vivo conditions than to the existing replica cast geometries. Furthermore, the numerical model of deposition appears to be capturing the correct relative magnitudes of molecular diffusion (proportional to d_p ) and convective diffusion (proportional to Q) based on comparisons to previous in vivo results.

Deposition by impaction is typically expected to correlate with the Stokes number, which captures the deviation of particles from curved streamlines. The relationship between inertial impaction and deposition has been clearly demonstrated for curved tubes (Cheng and Wang 1981), oral airways (Cheng et al. 1999; Xi and Longest 2007a), and respiratory bifurcations (Kim and Fisher 1994). The deposition of micrometer aerosols in the nasal geometry during inhalation has also been correlated with the Stokes number (Cheng 2003). However, Figure 9 clearly shows that the Stokes number alone (which is proportional to d ² _p Q) cannot be used to predict inertial deposition of particles ranging from 100 nm–1 μm in the nasal cavity during exhalation. Furthermore, this effect is observed in the laminar, transitional, and turbulent flow regimes. As a result, turbulent dispersions cannot fully account for the poor correspondence between deposition and St _k . Based on the numerical results of this study, the best correlation for inertial deposition was found to be using the impaction parameter d ^0.2 _p Q, which is proportional to a Stokes and Reynolds number combination of St ^0.1 _k Re ^0.9. As a result, the Stokes number appears to have a relatively minor role in the deposition of fine aerosols in the nasal cavity with the Reynolds number having a more dominate influence. This effect may be the result of significant changes in the flow field and strength of secondary flows as the flow rate is increased in the highly complex nasal geometry. In this case, the generation of new secondary flow features or a nonlinear increase in secondary flow magnitude contributes to deposition at a rate greater than the Stokes number effect for submicrometer aerosols. Similar results were observed by Longest and Vinchurkar (2009) for exhalation in a double bifurcation model.

Comparisons between the continuous-field drift flux approach and standard Lagrangian particle tracking indicate some differences in deposition for both ultrafine and fine respiratory aerosols. Considering the diffusion regime, predictions of the DF-VC model are consistently below the Lagrangian tracking results (Figure 6). These differences are a maximum for the smallest particles considered (approximately 8% difference) and decrease as the aerosol size increases. Longest and Xi (2007b) reported that non-continuum effects can influence the deposition of small nanoparticles in the range of 1–5 nm through two mechanisms. First, the inclusion of the Cunningham correction factor in the diffusion coefficient serves to increase diffusional effects. In addition, non-continuum flow effects can increase residence times of aerosols by increasing the slip velocity. The Lagrangian tracking model captures this effect in the drag term, whereas the DF-VC model approximates this phenomenon in the slip velocity calculation. It is expected that the increased residence times of small particles is better captured by the Lagrangian model resulting in an increase in deposition. Furthermore, the Lagrangian model is observed to better match the experimental data in the range of 1–5 nm (Figure 6). However, other factors in addition to increased residence times may also contribute to the higher deposition of 1–5 nm aerosols observed in the experiments (Cohen et al. 1995; Kelly et al. 2004b; Wang and Lai 2006). As a result, it is not clear if the Lagrangian or DF-VC model is more accurate for the prediction of ultrafine aerosol deposition.

Considering inertial effects, a comparison of the Lagrangian tracking approach and the DF-VC model indicates some differences in deposition (Figure 8). In each case, the Lagrangian model predicts a higher rate of deposition than with the DF-VC method. Both the Lagrangian and drift flux models predict significant increases in deposition compared to the CS model beyond St_k = 1 × 10⁻⁵ as a result of including particle inertial effects. However, it is not clear if the Lagrangian or DF-VC model better predicts aerosol deposition. Previous results indicate that the Lagrangian model typically over predicts inertial deposition due to turbulent dispersion when used with an isotropic turbulence model as employed in this study (Matida et al. 2004; Xi and Longest 2007a). Furthermore, Longest and Xi (2007b) showed that inadequacies in accounting for near-wall effects can lead to artificially elevated Lagrangian particle deposition. In contrast, the DF-VC model was previously validated to accurately predict the deposition of submicrometer aerosols in both the diffusional and inertial regimes. Furthermore, the DF-VC model was found to accurately predict millimeter-scale localized particle deposition and hotspot formation for submicrometer aerosols (Longest and Oldham 2008). In contrast, comparison of Lagrangian model results with a highly localized in vitro data set indicated significant differences in deposition patterns (Longest and Oldham 2006).

In this study, deposition results were presented primarily for the drift flux model based on previous validations of regional and highly localized deposition in a respiratory bifurcation model (Longest and Oldham 2008). However, the DF-VC model may under predict the deposition of submicrometer aerosols, as observed by Xi and Longest (2008c; 2008b). An expression comparable to Equation (14) for extrathoracic deposition during nasal exhalation based on Lagrangian model results can be written

Considering the observations of Matida et al. (2004) and Longest and Xi (2007b, 2007a) as described above, deposition predictions from this Lagrangian-based correlation may be higher than experimental observations. Further studies are needed to better validate both the Lagrangian and DF-VC models for determining aerosol deposition in complex respiratory geometries like the nasal cavity, especially for particles in the range of 100 nm and larger. However, the current drift flux (Equation [14]) and Lagrangian (Equation [16]) model correlations provide a reasonable range of potential deposition values.

One uncertainty that may produce differences between simulation results and in vivo data is the actual uvula position during normal expiratory breathing. The current pharyngeal–laryngeal model was based on CT images acquired with the patient in the supine position and at the end of exhalation. It is likely that the actual upper airway geometry over the exhalation period with an upright position is different from the current model. In particular, the uvula is a highly dynamic structure, with its position substantially determining the lower nasopharynx lumen, and thereby the inlet flow conditions into the nasal cavity during exhalation. In this study, the uvula is positioned on the back of the throat and obstructs a portion of the lower nasopharynx, which divides the lower nasopharynx into thin channels. As a result, elevated depositions were observed around the uvula-formed nasopharynx division which may also influence deposition in the downstream nasal cavity during exhalation.

As described, the flow and particle profiles entering the model at the larynx were assumed to be uniform based on dispersive mixing in the trachea. This approximation likely has the largest impact on the local deposition characteristics in the larynx. However, the laryngeal jet that forms during exhalation is expected to control both the flow field and particle profiles entering the remainder of the extrathoracic region. As a result, inclusion of more realistic downstream conditions is not expected to have a significant impact on local and region deposition in the structures upstream of the larynx during exhalation.

Factors that limit the physiological realism of the current study include the assumptions of steady flow, simplified inlet conditions, a smooth and rigid airway surface, no humidity, and a constant nasal valve for various breathing conditions. Other studies have highlighted the physical significance of transient breathing (Shi et al. 2006), inlet velocity profiles (Keyhani et al. 1995; Subramaniam et al. 1998), nasal wall motion (Fodil et al. 2005), nasal valve change during respiratory maneuvers (Bridger 1970; Bridger and Proctor 1970), and intersubject variability (Storey-Bishoff et al. 2008). Moreover, due to limited access to high-quality medical images that span the extrathoracic airways, the model in this study was based on images of two patients, which may be different from an airway geometry acquired from a single subject. Each of these factors influences the realism of the model predictions in relation to actual particle deposition in the nose. However, results of the scan-based geometry and numerical model considered in this study matched the available in vitro and in vivo submicrometer deposition data for nasal exhalation to a high degree.

In conclusion, submicrometer aerosol deposition has been evaluated in a scan-based nasal–laryngeal model during expiration using both DF-VC and Lagrangian methods. Simulation results provide a reasonably good match to available in vitro and in vivo experimental deposition data and corroborate the use of an existing in vivo–based deposition parameter proposed by Cheng (2003), with the extension of including the Cunningham correction factor. Based on a broad range of flow rates and particle sizes, a new diffusion correlation was developed for ultrafine particles. The diffusion correlation was further extended to include the inertia-associated effects of fine particles. Specifically, a correction was proposed for extrathoracic deposition during nasal exhalation that is applicable for both ultrafine and fine respiratory aerosols. Complementary in vivo and in vitro studies on fine particle deposition are needed to further evaluate these numerical predictions, thoroughly test the proposed deposition correlations, and better understand deposition mechanisms during exhalation. Further studies are needed to address tidal breathing, surface compliance, vocal cord movement, upstream flow effects, and inter-subject variability before current numerical nasal airway models can be directly applied to make dose–response and health effects predictions.

The authors thank Dr. Julia Kimbell of the Hamner Institutes for Health Sciences for providing the MRI tracing data used in constructing the nasal cavity model. The use of CT data provided by Dr. Karen A. Kurdziel in the VCU Department of Radiology and Molecular Imaging Center under IRB approval (06263) is gratefully acknowledged. This work was sponsored by Philip Morris USA.

Notes

∗ANOT: Adult-Nasal-Oral-Trachea model.