Particle deposition in the pulmonary region of the human lung: A semi-empirical model of single breath transport and deposition

Aerosol Science 38 (2007) 228 – 245 www.elsevier.com/locate/jaerosci

Particle deposition in the pulmonary region of the human lung: A semi-empirical model of single breath transport and deposition S.S. Parka,∗ , A.S. Wexlera, b, c a Department of Civil and Environmental Engineering, University of California, Davis, CA 95616, USA b Department of Mechanical and Aeronautical Engineering, University of California, Davis, CA 95616, USA c Department of Land, Air, and Water Resources, University of California, Davis, CA 95616, USA

Received 8 August 2006; received in revised form 16 November 2006; accepted 20 November 2006

Abstract Particle deposition has been observed in the distal pulmonary regions of human airways, where their clearance is lower and thus their potential health effects greater. Whereas gases are transported convectively to proximal pulmonary regions and then diffuse to the distal ones, particles have a very low diffusivity so only few can be transported to the two or three most distal generations by convection in a single tidal breath. Yet models of particle deposition in human airways invariably employ a single breath to characterize regional pulmonary deposition. To study particle transport and deposition during multiple breathing cycles, it is essential to develop a single breath model that describes aerosol bolus dispersion and deposition for the entire lung. Here we present a simple semi-empirical model for whole lung aerosol bolus dispersion, and use bolus dispersion data from the literature to identify model parameters. The model assumes that the particle transport on inhalation and exhalation differs due to secondary flow asymmetries at the bifurcations. By adjusting model parameters to measurements, the effective particle transport profiles averaged over the entire lung and mixing intensity at different lung depths are obtained. The average inspiratory particle transport profile is found to be slightly sharper than parabolic, whereas on exhalation it is relatively blunt. The mixing intensity due to secondary flows is large in the conducting airways (∼ 95 ml of lung depth) but drops sharply due to the rapid decrease in Dean number in the deeper lung depths. Therefore, significant mixing due to secondary flows occurs mainly in the proximal airways during exhalation. A set of parameters is identified such that modeled bolus standard deviation, mode shift, and skewness agree with the measurements. Particle loss due to deposition was calculated in each generation and the resulting predictions of particle deposition fractions were in a good agreement with measurements. Although this model does not simulate detailed flow mechanics in the lungs, it does represent overall particle transport phenomena that lead to bolus dispersion observations and it provides a framework for the simulation of pulmonary particle transport during multiple breaths, as presented in an accompanying study of this work. 䉷 2006 Elsevier Ltd. All rights reserved. Keywords: Aerosol bolus dispersion; Particle transport profile; Mixing intensity

1. Introduction Transport of gases into and out of the lung is caused by a combination of convection, primarily in the proximal airways, and diffusion, primarily in the distal airways (Heyder, Blanchard, Feldman, & Brain, 1988; Ultman, 1985; ∗ Corresponding author. Tel.: +1 530 754 4963; fax: +1 530 754 4962.

E-mail address: [email protected] (S.S. Park). 0021-8502/$ - see front matter 䉷 2006 Elsevier Ltd. All rights reserved. doi:10.1016/j.jaerosci.2006.11.009

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

229

Nomenclature C CMIi Dei di ex in Li LEVi LMIi MS n Ni Q r R Rc,i Rei RN S Skew STD t TV EV ui  v(r) v Vp V (r) z  (nin , nex , S)

aerosol concentration volumetric cumulative mixing intensity to generation i dean number in generation i diameter of the airway in generation i expiratory phase inspiratory phase length of the airway in generation i local effective volume in generation i volumetric local mixing intensity in generation i mode shift of the bolus exponent of effective particle transport profile equation number of airways in generation i flow rate radial coordinate in local lumen radius of the local lumen radius of the airway curvature in generation i Reynolds number in generation i random number generated from a standard normal distribution scaling factor in the equation for CMIi skewness of the bolus standard deviation of the bolus time tidal volume expired volume local flow axial velocity in generation i kinematic fluid viscosity axial velocity as a function of r average axial velocity penetration volume volumetric lung depth of particle as a function of r axial distance of particle ratio of r and R value of the sum of squares between measurements and predicted parameters as functions of nin , nex , and S

Ultman & Thomas, 1979). Since these processes occur simultaneously, experiments with gases do not elucidate the separate role of these two mechanisms. Particle transport in airways is mainly governed by convection since particles, especially those in the 0.1–1 m size range, have mobility orders of magnitude lower than gases (Darquenne & Prisk, 2003; Heyder et al., 1988; Schulz, Schulz, & Heyder, 1996). Thus aerosol bolus dispersion provides a probe into the mechanisms of non-diffusive transport in the airways. During single breath aerosol bolus dispersion experiments, a small volume of particle-laden air is inhaled and the concentration of the particles as a function of volume is monitored during the subsequent exhalation. That the exhaled particles are dispersed over a much larger volume can only be explained by mixing of the inhaled bolus with particle-free air retained in the airways from previous breaths. This leads to transport of inhaled particles into the deeper lung, thus having a profound influence on regional particle dosimetry (Brand, Rieger, Schulz, Beinert, & Heyder, 1997; Davies, Heyder, & Subba Ramu, 1972; Heyder et al., 1988; Sarangapani & Wexler, 1999, 2000; Taulbee, Yu, & Heyder, 1978). Some investigators have used an effective diffusivity (Da ) derived empirically from gas bolus experiments in fully developed laminar flow in a straight tube to predict particle transport in the lung (Darquenne, Brand, Heyder, & Paiva, 1997; Darquenne & Paiva, 1994; Taulbee et al., 1978; Ultman, 1985; Ultman & Thomas, 1979). In diffusion, the axial

230

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

variance (2 ≡ 2Da t) of a bolus increases linearly with processing time (t). Gas bolus dispersion experiments exhibit a linear relationship between volumetric variance and respired flow rate (Ultman & Thomas, 1979). However, particle bolus dispersion experiments do not demonstrate a linear variance–time relationship; many studies show that measured aerosol bolus dispersion in terms of standard deviation, , is linearly dependent upon the bolus penetration lung depth (Brown, Gerrity, Bennett, Kim, & House, 1995; Heyder et al., 1988; Rosenthal, Blanchard, & Anderson, 1992), so variance, 2 , for particles as a function of penetration volume is faster than linear (Henry, Butler, & Tsuda, 2002). Taylor dispersion leads to radial concentration gradients but this will not lead to radial particle transport since particle diffusivity is orders of magnitude lower than that of gases. Therefore, axial particle bolus spreading does not follow the basic rules described in classical gas diffusive transport theories. In fact, time appears to play a less significant role in particle dispersion where the variance is more typically related to inhaled volume, implying that the aerosol bolus dispersion processes are more related to the volume traversed than the elapsed time (Heyder et al., 1988). Due to intricate flow patterns and complexity of the human lung geometry, computational fluid dynamics (CFD) simulations for particle transport in the tracheobronchial tree have been limited to the study of simple bifurcations or only several generations (Comer, Kleinstreuer, & Kim, 2001; Comer, Kleinstreuer, & Zhang, 2001; Zhang & Kleinstreuer, 2002). Although CFD models provide exquisite detail of air flow patterns and particle transport, they are too computationally intensive for simulating many generations, length scales, breathing cycles, and the overall measured dispersion. In this study, we propose a semi-empirical model for aerosol bolus dispersion; free parameters in the model are adjusted to match the dispersion measurements of Brand et al. (1997). Scherer and Haselton (1982) showed that the difference in the shape of inspiratory and expiratory velocity profiles can explain bolus dispersion. This difference implies that fluid particles carried forward during inspiration will not return to their starting points at the end of expiration, so that over a flow cycle, particles will be displaced upstream and downstream from their starting positions. At bifurcations, flow visualization studies show minimal secondary motions on inhalation with distinct velocity peak near the inner wall as flow enters the daughter branches, while substantial secondary motions are observed on exhalation forming quadruple vortices as the flow combines in the parent airway even for low Reynolds number region (Schroter & Sudlow, 1969). Due to the relative lack of secondary motions during inhalation, particles are thought to stream into the airways (Briant & Lippmann, 1993; Scherer & Haselton, 1982; Zhao, Brunskill, & Lieber, 1997). On exhalation, the story is different where at each bifurcation, secondary flows affect each particle’s radial position. This change in particle radial position at each bifurcation causes it to experience a range of axial velocities because the axial velocity is a function of radial position. Therefore, in the upper airways, secondary swirling of the expiratory flow over a series of bifurcations makes velocity profiles flatter than inspiratory profiles. In the respiratory zone where Reynolds number falls below 1, velocity profiles are dominated by viscous forces. However, the acinar flow is not kinematically reversible (i.e., during exhalation each fluid particle does not retrace the path taken during inhalation) since rhythmically expanding and contracting alveolar wall motions cause a chaotic acinar flow resulting in substantial flow irreversibility (Henry et al., 2002). Hence, the axial transport of a particle on exhalation in the pulmonary region might not correspond to low-Reynolds number laminar flow profile. In this model, we assume that the primary causes of bolus dispersion are asymmetries between the effective particle transport profiles during inhalation and exhalation, as proposed by Scherer and Haselton (1982), combined with secondary motion intensity as a function of lung depth characterized by local geometry and flow characteristics. These mechanisms are incorporated into a model that simulates both particle transport by convective flow in the lumen and mixing at bifurcations to characterize single breath aerosol bolus dispersion in terms of standard deviation, mode shift (MS), and skewness, which are commonly measured descriptors of width and shape of an aerosol bolus (Brand et al., 1997; Brown et al., 1995; Darquenne et al., 1997; Darquenne & Paiva, 1994; Scheuch & Stahlhofen, 1994). The obtained effective particle transport profiles on inhalation and exhalation represent an averaged axial particle transport over the whole lung during inhalation and exhalation, respectively. 2. Model development 2.1. Asymmetry in the particle transport profiles between inhalation and exhalation Consider symmetric well-developed laminar flow profiles during inhalation and exhalation. A thin aerosol bolus is stretched into parabolic surface of rotation during inhalation (Henderson, Chillingworth, & Whitney, 1915). Since the flow is reversible, the stretched bolus will reassemble during exhalation resulting in zero dispersion as long as the

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

231

particle’s intrinsic motion is negligible. Hence, the difference in flow profile between inhalation and exhalation has been proposed as the primary cause for aerosol bolus dispersion (Briant & Lippmann, 1993; Sarangapani & Wexler, 1999; Scherer & Haselton, 1982). In this work, it is also assumed that this flow irreversibility causes dispersion in the human airways. This assumption is primarily a result of the nature of the secondary fluid flow in a bifurcating geometry (Sarangapani & Wexler, 1999; Schroter & Sudlow, 1969). Sarangapani and Wexler (1999) simulated secondary flow patterns in a bifurcating geometry at three different low Re values (1, 10, and 100) that are more appropriate to the conducting and pulmonary airways by creating a three-dimensional bifurcation model using FIDAP (Fluid Dynamics Analysis Package). A weak double vortex near the bifurcation was observed during inhalation and a stronger quadruple vortex was seen during exhalation. They calculated the ratio of radial to axial velocity close to the bifurcation to quantify the intensity of radial mixing. During inhalation, the ratio was less than 1% for flow simulations at Re = 1, slightly greater than 1% for Re = 10, and averaged at about 3% for flow simulation conducted at Re = 100. Based on these flow simulations, radial mixing due to secondary motion during inhalation in the vast majority of the human lung is insignificant compared to the influence of secondary motions during exhalation. Such small secondary motions during inhalation result in axial streaming of the inhaled particles into the daughter branches, causing inhaled particles to penetrate to deeper regions of the lung than would be expected by the mean fluid front (Briant & Lippmann, 1993; Scherer & Haselton, 1982), whereas on exhalation, the effective velocity profile of the particles is blunter due to secondary motions at the bifurcations that cause particles to change their radial position in the lumen as they pass each bifurcation (Briant & Lippmann, 1993; Fresconi, Wexler, & Prasad, 2003; Haselton & Scherer, 1982; Schroter & Sudlow, 1969; Zhao et al., 1997). In the pulmonary regions, even though the Reynolds number is low, alveolar flow can be chaotic due to rhythmical alveolar wall motions. Henry et al. (2002) showed that the shape of a tracer bolus evolves to a stretch-and-fold fractal-like pattern suggesting that aerosol transport in the alveolated region is primarily due to chaotic acinar flow irreversibility. Therefore, we assume that the irreversibility of the particle transport profile is caused by an asymmetry in the intensity of the secondary motions at each bifurcation during inhalation and exhalation. In this work and that of the measurements that we rely on, ∼ 0.5 m particles are considered due to their low diffusivity, so that they do not diffuse significantly during a tide, and their low Stokes number, so that their inertia does not lead to significant departure from flow streamlines. Thus the particles may be considered so-called fluid particles that track the flow. Although the Reynolds number is very low in most of the lung, the local particle transport profile is not necessarily parabolic due to a number of mechanisms that disrupt the flow. For instance, the centerline maximum repeatedly strikes carinas, slowing them down, and leading to blunter than parabolic flow profiles on average. Sharper velocity profiles in the upper airways are observed in association with the laryngeal jet and during the transition region after entering daughter airways (Simone & Ultman, 1982). In well-developed laminar flow in a lumen of circular cross section, the flow profile is parabolic so the ratio of maximum velocity at the center to the lumen averaged velocity is 2. Due to glottal constriction, the axial velocity is greatly accelerated so that the ratio of maximum velocity to mean velocity becomes 3.5. After this jet like behavior, about 10 cm downstream of the glottis, this sharp velocity profile experiences a rapid transition to a blunter profile where the ratio of centerline to average velocity is near 1.4, due to turbulent flattening of velocity profile (Schulz et al., 1996). On inhalation, the divergent flow profile in the bronchial tree could make the profile sharper than parabolic (Briant & Lippmann, 1993), but encounters with the carina are expected to blunt the profile. In the pulmonary regions, the alveolated structure of acinar ducts and rhythmical motion of the alveolar walls are likely to play a significant role in particle transport. The curvature of flow streamlines due to the presence of side-wall alveoli may cause the particles to enter the alveoli and deposit there. Recently, Haber, Yitzhak, and Tsuda (2003) investigated the effect of the rhythmical alveolar wall motion on the flow. They showed highly complex particle trajectories in the acinus under low Reynolds number flows in a cyclically expanding and contracting alveolated duct structure (Haber et al., 2003). The particle transport profile does not represent these flow details but instead represents the transport averaged over entire lung. Aerosol bolus experiments are non-invasive so only characterize the shape of the bolus averaged over the lumen before inhalation and after exhalation. The exhaled particle concentration profile is determined by all the transport processes integrated over the inhalation and exhalation tidal cycle. We assume that inspiratory and expiratory effective particle transport profiles can be used to characterize the particle transport velocities on inhalation and exhalation, respectively. Steady particle transport profiles are assumed because during normal breathing (flow at the mouth of 500 ml/s, at the frequency of 0.25 Hz), the Reynolds number is less than approximately 400 and the Womersley number () is less than about 0.6 in the regions below of the fourth generation. Since  and flow amplitude decreased rapidly after the first seven generations, quasi-steady flow can be assumed in most regions of the lung (Grotberg, 2001;

232

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

A(z) (z) Axial distance, z Volumetric lung dept, Vinproximal Volumetric lung depth,Vindistal Fig. 1. Trumpet model of airways.

Ramuzat & Riethmuller, 2002). In summary, the flow profile in the airways is typically laminar, developing near the bifurcations, and has small secondary motions that are stronger on exhalation than inhalation. Bolus dispersion measurements characterize particle transport averaged over the entire lung and this transport can be characterized by effective particle transport profiles that differ between inhalation and exhalation. Due to secondary motions at the bifurcations in the upper airways and chaotic acinar flow in the pulmonary regions, the particles change their radial position so sample different axial velocities as they move from generation to generation. Thus the particle velocity and effective transport profiles may not correspond, but it is the transport profiles that govern dispersion in the airways, so from this point onwards will confine our discussion to transport profiles. The local transport properties in an airway are duplicated in numerous parallel airways at the same generation, so we characterize the a simple effective particle transport profile in a trumpet (Yu, 1978) whose bluntness depends on the variable n    r n    r n   n + 2  Q(t)  n+2 v(r, z, t) = = , (1) · v(z, t) · 1 − · · 1− n R n A(z) R where v is the local axial transport of particles as a function of radial position, r, the local lumen radius is R, v is the local average velocity defined as the volume flow rate, Q(t), divided by its cross sectional area, A(z), at the lung depth z, as illustrated in Fig. 1 (Weibel, 1963). This profile satisfies the zero velocity boundary condition at the lumen wall and conservation of mass within the lumen in that the integral of the velocity over the lumen cross section gives the volumetric flow rate. In Eq. (1), n describes the bluntness of the profile with n = 2 corresponding to the parabolic profile characteristic of well developed laminar flow, while n > 2 characterizes average transport profiles that are blunter and n < 2 describes sharper transport profiles. The volumetric lung depth, V (z, r), at axial position z, is obtained by integrating A(z) over z:  V (z, r) =

z

A(z ) dz =

0



t 0

A(z )

dz  dt = dt 



t

A(z )v(r) dt  ,

(2)

0

where the volumetric lung depth depends on effective radial position, r, because the transport profile depends on radius. By substituting for v(r) from Eq. (1), A(z) cancels and Q(t) is the only quantity that is a function of time in Eq. (2). Therefore, this integral over inhaled time gives Eq. (3), the volumetric lung depth of the aerosol bolus as a function of radial position after inhalation, proximal

Vin

(r) =

  r nin  nin + 2 ·TV · 1− nin R

(3)

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

233

Particle-free area Particle-laden area

Exit, Vp=0

Fig. 2. Cartoon of an aerosol bolus after passing through multi bifurcations.

and during exhalation Vex (r) = Vin (r) −

  r nex  nex + 2 , · EV · 1 − R nex

(4)

where TV is the inhaled (tidal) volume, EV = Qt is the exhaled volume where t is time during the exhalation portion of the tidal cycle, and nin and nex are exponents describing the average effective inspiratory and expiratory particle transport profiles, respectively. To simulate aerosol bolus dispersion, the proximal and distal edges of the bolus are tracked as shown in Fig. 1. In the experiments of Brand et al. (1997), the particles are uniformly distributed within a narrow 20 ml bolus volume upon inhalation, so the proximal and distal edges of the bolus are separated by this volume in our bolus simulations. These edges are severely distorted as the bolus is inhaled and subsequently exhaled as shown in Fig. 2. Although the actual bolus shape at the exit may be quite contorted, the boundaries are still distinct because ∼ 0.5 m particles have very low diffusivity and therefore track fluid parcels as they are inhaled and exhaled. Hence at the mouth (Vp = 0) during exhalation, particles are only encountered between these proximal and distal edges, and if there was no loss due to deposition, the particle concentration within this volume would be equal to the initial concentration. During inhalation, particle-laden air streams into the airways with an effective transport profile exponent of nin . The particle deposition fraction in each generation on inhalation is calculated as described in the next section. After inhalation, the volumetric position of the distal edge (Vindistal ) is acquired by adding the bolus volume, 20 ml to the volumetric position proximal proximal ) in Eq. (3). Vindistal and Vin then describe the position of the edges of the bolus at of the proximal edge (Vin the beginning of exhalation. As particles are transported by an effective transport profile exponent of nex , deposition particle loss on exhalation is calculated in each generation. During typical bolus dispersion experiments, the subject’s standard tidal volume is inhaled but the subject fully exhales so that the full shape of the bolus can be measured. In Eq. (3), the inhaled volume is expressed in terms of tidal volume and radial position to characterize the position of the bolus at the end of inhalation, while Eq. (4) is expressed in terms of expired volume (EV) so that this full bolus shape during exhalation can be described.

234

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

2.2. Deposition mechanisms Particle deposition depends on morphometric and flow parameters. In this work, we develop a generation-bygeneration compartment model of the lung to estimate the particle deposition in the each generational compartment as the air traverses the airway network. Oral deposition efficiency as a function of Stokes number for an adult airway replica suggested by Cheng (2003) was applied. Equal deposition efficiency for oral deposition on inhalation and exhalation was assumed due to the absence of specific data (Cheng, Zhou, & Chen, 1999). Our thoracic lung compartment model consists of 24 generations, the first 15 generations forming the conducting airways and the rest the intra-acinar region. In the intra-acinar region, about 25% of the volume is contained in the lumen and 75% in the alveoli (Haefeli-Bleuer & Weibel, 1988; Weibel, 1963). Particle loss in each compartment of the tracheobronchial and pulmonary region is due to impaction, diffusion, and sedimentation. Empirical deposition equations for impaction (Kim, Iglesias, & Garcia, 1989a; Kim, Iglesias, & Garcia, 1989b), diffusion (Ingham, 1975) and sedimentation (Pich, 1972) are applied at each generational compartment to compute deposition. Inhaled particles in the range of 0.1–1 m are likely to be transported to the distal pulmonary regions of the lung since these particles are too large to undergo significant diffusion and too light to be substantially affected by inertia (Davies et al., 1972; Schulz et al., 1996). Many studies have explored gravitational settling and deposition of fine particles in the pulmonary regions, but most of them ignored or simplified the complexity of the flow that exists inside the alveoli. Haber et al. (2003) showed that alveolar wall motion is crucial in determining the enhancement of aerosol deposition inside the alveoli. Gravitational deposition processes and deposition sites in rhythmically expanding and contracting alveoli differ from conventional predictions made by classical models, which treat the acinar duct as a straight pipe with rigid walls. They suggested local alveolar deposition probability of the rhythmically moving and rigid-walled alveoli due to gravitational deposition for 0.5 m particles in generation 16–23. These gravitational deposition efficiencies in the pulmonary region are used instead of classical sedimentation efficiencies (Pich, 1972) to include the effect of complex acinar flow to the enhanced pulmonary deposition. Moving-wall deposition efficiency by Haber et al. (2003) might be an overestimate since deposition efficiencies for the moving alveoli wall are only relevant for horizontally placed alveolus with its mouth facing upward. Considering that about one quarter of the alveoli face upwards, we use only one fourth of Haber’s moving-wall deposition efficiencies in generation 16–23. Intraacinar airways are of three types: respiratory bronchioles, where only part of the surface is covered by alveoli; alveolar ducts, which are completely alveolated; and alveolar sacs, the terminations of alveolar ducts (Haefeli-Bleuer & Weibel, 1988). Since the transition to alveolar ducts is gradual, the first few generations of respiratory bronchioles (generation 14–16) are partly alveolated. In order to apply Haber’s mechanism for gravitational sedimentation deposition efficiency, measurement data of the human respiratory acinus, such as alveoli per duct, number of ducts, diameter and length of ducts, and diameter of alveoli in generation 14–23, are used to calculate the surface area ratio of alveolated portion to the non-alveolated region in each generation (Schreider & Raabe, 1981). By employing the concept of lung compartmental model, we track particle deposition in a generation as the particleladen air traverses each compartment. Particle loss in each compartment was calculated due to the combined action of all the deposition mechanisms. First of all, particle deposition in the upper respiratory tract was calculated by Eq. (5). C1 = C0 (1 − 0 ),

(5)

where C0 is the initial uniform particle concentration, C1 is the uniform particle concentration at the first compartment entrance, and 0 is oral deposition efficiency suggested by Cheng (2003). During inhalation, particles are transported by an effective transport profile exponent of nin . Depending on the locations of particles at the end of inhalation, particle deposition and concentration after particle deposition at each compartment were calculated by Ci (Vp , ) = C1

i−1

j

j

(1 − SD )(1 − I ),

(6)

j =1

where Ci is the uniform particle concentration at radial position, , and axial position, Vp in the ith compartment, SD is combined particle deposition efficiency due to sedimentation and diffusion, and I is particle deposition due to impaction. During exhalation, the position of each particle is tracked. Consider a particle that is located at radial position, , and axial position, Vp in generation 16 after inhalation. If this particle is transported to Vp in generation 14 after exhaling 1 ml, the concentration at radial position, , axial position, Vp in generation 14 is calculated

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

235

using C14 (Vp , ) = C16 (Vp , )

16

j

j

(1 − SD )(1 − I ).

(7)

j =15

The deposition fraction at a given volumetric lung depth (from 20 to 800 ml) is obtained as a sum of particle loss during inhalation and exhalation. Secondary motions at the bifurcations and chaotic acinar flow cause fluid particles to be transported to different radial positions and these different radial positions have different axial velocities from the ones expected by averaged flow profiles. As fluid particles traverse numerous bifurcations on inhalation or exhalation, the particle axial displacements have some variance about the means described by Eqs. (3) and (4) due to chaotic advection caused by secondary motions and alveolar wall motions (Haber et al., 2003). Modeling this variance about the mean position is discussed next. 2.3. Variance about the mean volumetric position due to radial motions During expiration, flow in the two daughter branches combine at each bifurcation to form a quadruple vortex pattern in the upper airways (Haselton & Scherer, 1982; Schroter & Sudlow, 1969; Zhao et al., 1997). Such motions lead to changes in radial and azimuthal positions of particles, which lead to changes in their axial transport. In the acinar flow, a tracer bolus undergoes cyclic stretch-and-fold deformation leading to chaotic mixing (Henry et al., 2002). If we consider a particle bolus that must pass many generations during exhalation, when it reaches the mouth, some particles will be ahead and some behind due to the numerous eddy motions, and consequent axial velocities, that they have experienced. Thus, although Eq. (4) expresses the overall average profile due to these motions on exhalation, there will be some chaotic variance about this average such that at a given exhalation volume, EV, some radial positions may have particle laden air while others will have particle free air at the exit, as shown in Fig. 2. As a descriptor of bolus dispersion, the standard deviation, which represents the width of expired particle concentration distribution, results from both the interaction between the average inhaled and exhaled transport profiles and the variance about these averaged positions due to secondary motions. Although we have represented the average position of the fluid particles by Eqs. (3) and (4), axial displacements from the mean positions of particles due to secondary motions in the upper airways and chaotic acinar flow in the pulmonary regions need to be described to simulate aerosol bolus measurements. The intensity of secondary motion at each volumetric lung depth is modeled by taking into account the concept of local effective volume based on Weibel’s model (Weibel, 1963) and the local flow properties characterized by the local Dean number. In the airways, individual branch diameters decrease with increasing generation, while the corresponding values of summed cross-sectional area increase (Hammersley & Olson, 1992; Horsfield, Dart, Olson, Filley, & Cumming, 1971; Ultman, 1985; Weibel, 1963). In particular, the total cross-sectional area of the airways in the lung periphery (i.e., about 200 ml in volumetric lung depth) rapidly increases leading to a decrease in axial gas velocity. Even though flow velocity decreases dramatically beyond the conducting airways and is insignificant in the peripheral lung, it has been shown that convective mixing, which refers to bulk processes that irreversibly transfer particles from the tidal air to residual air, occurs in the periphery in similar magnitude to that in the central airways from nonintrinsic particle bolus dispersion experiments (Heyder et al., 1988; Taulbee & Yu, 1975; Taulbee et al., 1978). Therefore, lung geometrical factors such as the large number of airway branches at higher generation may influence the intensity of secondary motions in these deeper regions. Although the intensity of secondary motions at one bifurcation is weak in the low Reynolds number region, particles experience many bifurcations as they flow into and out of the pulmonary region. During exhalation, particles in the more distal regions pass through more bifurcations before they reach the mouth so it is expected that the intensity of secondary motions will be proportional to the integrated effect of all bifurcations that each particle passes during exhalation from its position after inhalation back to the mouth. We assumed that the level of secondary motion can be characterized by geometric factors and flow patterns at each volumetric lung depth. There are two factors that influence the magnitude of secondary motion at a given volumetric lung depth: the local effective volume and the Dean number dependent mixing intensity. Let us first characterize the local effective volume (LEVi ) at generation i. The quadruple vortices formed at each bifurcation on exhalation have a size that is governed by the lumen diameter, di —larger diameter results in larger

236

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

vorticies. But during synchronous ventilation, this mixing occurs at many bifurcations simultaneously at the same volumetric depth. The local volume occupied by these vortices at generation i is the product of this characteristic vortex size, di , and the overall cross sectional area, A(zi ) = Ni di2 /4, giving    d 2 di = Ni di3 , LEVi = Ni (8) 4 i 4 where Ni is the number of branches and di is the average lumen diameter at generation i. On exhalation at each bifurcation, the eddies swirl in the local volume, LEVi , with an intensity that is proportional to the flow magnitude and airway dimensions via the Dean number, a dimensionless group characterizing secondary flow in curved channels, defined by       Li 0.5 Li 0.5 u i di Dei = = Rei · , (9) ·  2Rc,i 2Rc,i where ui is local flow velocity, Li is the airway length, Rc,i is the airway radius of curvature, at each generation i, and  is fluid kinematic viscosity. De is proportional to Re, which represents the ratio of inertial to viscous forces in the flow. Although the bifurcating structure of pulmonary airways is complex and there is a large variation in their geometric details, the ratio of airway length to radius of curvature in each generation is approximately 1.2 as employed in both asymmetrical and symmetrical branching lung models (Darquenne, 2001; Horsfield et al., 1971; Weibel, 1963). Therefore, it is assumed that the ratio of airway length to radius of curvature in each generation is roughly invariant with generation, leading to the conclusion that the magnitude of secondary motions at a given lung depth is proportional to the local Reynolds number. Thus we characterize the local mixing intensity (LMIi ) in terms of volume influenced by these two main factors—local effective volume (LEVi ) and Reynolds number (Rei )—given by LMIi = LEVi × Rei =

Q 2 d .  i

(10)

As is clear from Eq. (10), LMI drops rapidly with generation i due to decreasing di , so that the local mixing intensity will be insignificantly small in the pulmonary airways. During exhalation, particles experience secondary motions at each bifurcation. We assume that the transport effect of these motions is proportional to LMIi . Particles that are inhaled to deeper volumes pass more generations during exhalation so experience more secondary motions. For instance, the tip of the aerosol bolus, located at the deepest lung depth at the end of inhalation, experiences more bifurcations during exhalation and therefore more secondary swirling motions than other particle positions in the bolus. From their position at the end of inhalation, particles pass all the local mixing volumes on their way back to the mouth. Therefore, the variance in bolus position is proportional to the integral of Eq. (10) from its deepest volumetric lung depth (k = i) to the mouth (k = 1) during exhalation. CMIi = S · RN() ·

i

k=1

Q 2 dk ,  i

LMIk = S · RN() ·

(11)

k=1

where CMIi is the cumulated mixing intensity for parcel of air starting at generation i, S is a scaling factor that will be identified from the measurement data to acquire the range of mixing intensities at each lung depth and RN() is a random number as a function of normalized radial position to simulate the pattern of chaotic advection,  = r/R, generated from a standard normal distribution with a zero mean and a standard deviation of one to simulate axial variations on the averaged transport profile due to secondary motions. Therefore, CMIi is a function of radial and axial position. This axial mixing volume (CMIi ) is added to the mean position of each particle after inhalation. 3. Parameter identification and model evaluation The model just described was used to simulate the same processes as experiments of Brand et al. (1997), who employed 20 ml of aerosol bolus at the beginning of inhalation of ∼ 0.84 m diameter aerosol particles. The breathing condition is simulated at a constant flow rate of 250 ml/s and different mean penetration volumes, from 20 to 800 ml in

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

237

the aerosol bolus experiments. The bolus is stretched by inspiratory flow with a velocity profile characterized by nin . The penetration volume of the bolus is determined by Eq. (3). On exhalation, the cumulative mixing intensities (CMIi ) are applied depending on the volumetric positions of particles. Local mixing intensities (LMIi ) in Eq. (10) range from 4850 ml at the first generation to ∼ 2.5 ml at the most distal generations. Although the local mixing intensity is low in most airways, on exhalation particles starting at deeper regions pass through many bifurcations before they reach the mouth so it is the integrated effect of passing many bifurcations that is responsible for the overall mixing intensity. Due to the very small Dean number in the pulmonary regions, the local mixing effect is small. Particles are exhaled with an effective particle transport profile governed by the value of nex . The ratio of the observed concentration at the mouth at given time (t) to the particle concentration in the inhaled bolus, Ct, /C0, , is given by averaging over the lumen. Ct =

1 · C(Vp , ) · 



1

2 d =

1

0

2 () · C(Vp , ) ·  · ,

(12)

0

where  = 0.0125 is the interval of radius ratio. () is a function that is zero where there are no particles and one where there are particles at each normalized radial position, . C(Vp , ) is the particle concentration after deposition as a function of the particle’s axial and radial positions. Parameters such as standard deviation, half width, MS and skewness are quantified to measure the broadening of the aerosol bolus in the respired air in the same calculation process as Brand et al. (1997). Aerosol bolus descriptors were calculated in terms of standard deviation, MS, and skewness using the standard formulas as specified in Brand et al. (1997). Integration of the inhaled and exhaled bolus involved aerosol concentrations that are higher than 15% of the exhaled bolus peak concentration as recommended by Brand et al. (1997). In brief, standard deviation () is calculated by the square root of the difference of the second moment of exhaled and inhaled boluses (2ex − 2in ). MS is obtained by the difference between the volumetric lung depth (Vp ) of the inhaled bolus and bolus mode (M), which is the cumulative exhaled air volume at which the maximum particle concentration is measured. As for skewness, it is computed as the third moment of the exhaled bolus divided by 3 . All the bolus descriptors, such as standard deviation, MS, and skewness, are used to adjust model parameters (i.e., nin , nex , S) to the measurement data. The values of these free parameters are obtained by minimizing the sum of the differences between measurements and predicted results at 8 mean penetration volumes from 20 to 800 ml. Eq. (13) is described as follows:

(nin , nex , S) ⎡ 2 2 800

ml STDP ,Vp (nin , nex , S) MSP ,Vp (nin , nex , S) ⎣ 1− = + 1− STDM,Vp MSM,Vp Vp =20 ml

+ 1−

SkewP ,Vp (nin , nex , S) SkewM,Vp

2 ⎤ ⎦,

(13)

where STDM,Vp and STDP ,Vp are standard deviations at the given mean penetration volumes from measurements and simulations, respectively, and MS and Skew are the corresponding MS and skewness, respectively. The optimized values of these parameters (i.e., nin , nex , S) and comparison to aerosol bolus dispersion measurements are described in the next section. 4. Results Minimizing Eq. (13) yields nin = 1.7 and nex = 5.0 for the effective particle transport profile exponents and S = 0.01 for the scaling factor of cumulative mixing intensity. Fig. 3 illustrates that asymmetry in the particle transport profiles between inhalation and exhalation results in the observed aerosol bolus dispersion. If the inspiratory transport profile is the same as the expiratory one, it yields no dispersion (as shown in Fig. 3A), whereas particles disperse two times the mean penetration volume (500 ml) with a parabolic inspiratory transport profile and a plug expiratory transport profile (as shown in Fig. 3B). When the mixing intensity is superimposed on the different transport profiles during inhalation

238

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

Particle Concentration ratio

1.0

0.8

0.6

0.4 Aerosol bolus 20 ml

0.2

0.0 0

500 Inhaled volume (ml)

1000

1.0

0.10

0.08

0.8

0.06

0.6 (A)

Particle Concentration ratio

A

0.4

0.04 D 0.02

0.2

B C

0.0

0.00 0

500 Exhaled volume (ml)

1000

Fig. 3. Exhaled particle concentration distributions with different set of effective particle transport profiles; A: nin = 2 on inhalation, nex = 2 on exhalation, S = 0 (no mixing), B: nin = 2 on inhalation, nex = ∞ (plug flow) on exhalation, S = 0 (no mixing), C: nin = 2 on inhalation, nex = ∞ (plug flow) on exhalation, S = 0.01, D: nin = 1.7 on inhalation, nex = 5 on exhalation, S = 0.01.

and exhalation (S = 0.01), the exhaled particle concentration distribution is no longer smooth representing the pattern of the observations and particles are more dispersed compared to the results without mixing (as shown in Fig. 3C). The exhaled particle concentration with the optimized parameter values (nin = 1.7, nex = 5.0, S = 0.01) is shown in Fig. 3D. Exhaled particles are less dispersed compared to Fig. 3C. The effective inspiratory particle transport profile with the exponent value nin = 1.7 is near parabolic yet slightly sharper than laminar velocity profile indicating that during inhalation the secondary motions are weak, the Reynolds numbers are low in most of the lung volume, and the streaming is large—larger than for a parabolic profile. However, the effective expiratory particle transport profile is not plug as previous researchers assumed (Sarangapani & Wexler, 1999, 2000), but is more dispersive. When nex = 5.0, the tip of bolus penetrates 1.4 times the mean penetration volume during exhalation compared to the penetration volume when nex = ∞ corresponding to plug flow. Bolus description parameters (MS, STD) resulting from different sets of nin , nex , and S, as shown in Fig. 3, were compared to observations in Fig. 4. When the effective expiratory transport profile is equal to the inspiratory one, MS and STD is constant over all penetration lung depths since there is zero dispersion (Fig. 4A). With parabolic inspiratory

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

239

50 A, B Bolus mode shift (ml)

0

C

-50

-100

D

-150

-200 0

200

400 600 Volumetric lung depth (ml)

800

Bolus standard deviation (ml)

500

400

300

200

C B D

100

A 0 0

200

400 600 Volumetric lung depth (ml)

800

Fig. 4. Comparison of predicted parameters resulted from different set of nin , nex , and S to the measurements; A: nin = 2 on inhalation, nex = 2 on exhalation, S = 0 (no mixing), B: nin = 2 on inhalation, nex = ∞ (plug flow) on exhalation, S = 0 (no mixing), C: nin = 2 on inhalation, nex = ∞ (plug flow) on exhalation, S = 0.01, D: nin = 1.7 on inhalation, nex = 5 on exhalation, S = 0.01. • (mean ± SD) refers to measurement data from Brand et al. (1997).

transport profile and a plug expiratory one without mixing, particles dispersion increases with the penetration lung depth showing a linear increase in STD above the measurement data in Fig. 4B. Since the maximum concentrations can be considered as medians in the exhaled particle concentrations, MS as a function of Vp becomes constant. When mixing was applied to the case of Fig. 4B, particles disperse more resulting in a small increase in STD especially at lower Vp (< 300 ml) (Fig. 4C) while the range of MS is similar to Fig. 4B. With the optimized parameter values, both of the modeled STD and MS as a function of Vp are in good agreement with the measurements in Fig. 4D. Cumulative mixing intensities were calculated at each penetration volume and the scaling factor value, S = 0.01, was identified from the minimization of Eq. (13). CMIi as a function of Vp is shown in Fig. 5. Particle dislocations from the obtained transport profiles can be described in terms of volumetric lung depth. For example, a particle located at 600 ml after inhalation may experience a chaotic axial dislocation of about 95 ml of volumetric lung depth from the mean particle transport profile during exhalation. The cumulative mixing intensity given by Eq. (11) sharply increases in the conducting airways (∼ 95 ml) but levels off quickly due to the rapid decrease in Dean number with volumetric depth. Therefore, only the proximal airways contribute to the cumulative mixing intensity. Several researchers found that as the aerosol bolus penetrates deeper into the lung, it becomes more dispersed (Brand et al., 1997; Brown et al., 1995; Darquenne et al., 1997; Darquenne & Paiva, 1994; Heyder et al., 1988; Schulz et al.,

240

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

Cumulative mixing intensity (ml)

100 80 60 40 20 0 0

40

80 120 Volumetric lung depth (ml)

160

Fig. 5. Volumetric cumulative mixing intensity as a function of volumetric lung depth with optimized scaling factor, S = 0.01.

Aerosol concentration (C/Co)

0.12 Vp=100 ml

0.10

Vp=200 ml 0.08 Vp=400 ml 0.06

Vp=600 ml Vp=800 ml

0.04 0.02 0.00 0

500

1000 1500 Expired volume (ml)

2000

2500

Fig. 6. Particle concentration distributions as functions of the expired air volume and mean penetration volume (Vp ).

1994; Schulz et al., 1994). The concentration distributions predicted by this model are comparable to those measured by others (Darquenne et al., 1997; Heyder et al., 1988). The predicted concentration distributions of exhaled particles are shown in Fig. 6. The bolus becomes dispersed and the heights of maximum concentrations are reduced as bolus penetrates deeper regions in the lung. Due to the variance in axial velocity caused by secondary motions, the concentration distributions are not smooth. To reduce the effect of noise in the concentration distributions, the volumetric centre of mass for the upper 90% of the concentration distribution is calculated to identify the maximum concentration point for calculation of MS. Comparisons between the predicted and other measured bolus description parameters are displayed in Fig. 7 (Anderson, Hardy, Gann, Cole, & Hiller, 1994; Brand et al., 1997; Darquenne et al., 1997; Siekmeier, Schiller-Scotland, Kronenberger, & Gebhart, 1992; Verbanck, Schuermans, & Vincken, 2001). Since model parameters are optimized for satisfying three different bolus descriptors (Brand et al., 1997), the results of comparison does not fit well with each of the measured ones. In addition, there is a substantial intersubject variability in the data and the current model does not capture subject-specific parameter. Considering significant variations between measurements, which are even greater than the intersubject variability, our predictions are comparable and better fitted to the measurements by Brand et al. (1997) compared to other model prediction (Darquenne et al., 1997; Sarangapani & Wexler, 1999), and do show similar

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

A

B

50

300 Bolus standard deviation (ml)

Bolus mode shift (ml)

241

0 -50

-100

Our prediction Brand et al. (1997) Siekmeier et al. (1992)

-150

Anderson et al. (1989)

-200

250 200 150 Our prediction Brand et al. (1997)

100

Siekmeier et al. (1992) Schulz et al. (1995) Verbanck et al. (2001)

50 0

0

200 400 600 Volumetric lung depth (ml)

800

0

200 400 600 Volumetric lung depth (ml)

800

C

Bolus skewness

1.5

Our prediction Brand et al. (1997) Darquenne et al. (1997) Siekmeier et al. (1992) Schulz et al. (1995) Verbanck et al. (2001)

1.0

0.5

0.0

-0.5 0

400

800

1200

1600

Volumetric lung depth (ml)

Fig. 7. Predicted parameters optimized to measurement data, A: mode shift, B: standard deviation, C: skewness; • (mean ±SD) refers to measurement data from Brand et al. (1997), , predicted results from the model for nin =1.7, nex =5.0, S =0.01. Comparison of other experimental data (Anderson et al., 1994; Darquenne et al., 1997; Schulz et al., 1995; Siekmeier et al., 1992; Verbanck et al., 2001) is also demonstrated.

trends to all other available data (Anderson et al., 1994; Darquenne et al., 1997; Schulz et al., 1995; Siekmeier et al., 1992; Verbanck et al., 2001). The measured MS is close to zero for shallow lung depths, but decreases sharply as a function of Vp after 100 ml reaching −70 ml for Vp = 800 ml (Brand et al., 1997). The predicted MS agrees well with this trend considering the high intersubject variability (Fig. 7A). Fig. 7B shows the predicted and measured standard deviations. Both increase with increasing Vp showing nonlinear characteristics particularly for the shallow lung depths (Vp < 150 ml). The predictions agree with measurements of Scheuch and Stahlhofen (1991) showing that the increase in dispersion during shallow lung depths (20 < Vp < 100 ml) is steeper than in distal ones. Considering that skewness is a higher moment measure of bolus shape than standard deviation or MS, the predictions of skewness are agreed well with the observations, as shown in Fig. 7C. The predictions are slightly underestimated compared to measurements except for at the volumetric lung depth of 800 ml. Since deposition parameters are not adjusted by the measurements in this model, the model can be validated by comparing the predicted depositions to the measured ones. The predicted deposition fractions are compared to the experimental data in Fig. 8. The deposition fraction obtained by the sum of particle loss during inhalation and exhalation is shown as the ratio of deposited particles to total inhaled. They are in a good agreement considering the intersubject variability on the measurement data.

242

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

50

Deposition (%)

40

30

20

10

0 0

200

400 600 Volumetric lung depth (ml)

800

Fig. 8. Comparison of deposition fraction; • (mean ± SD) refers to measurement data from Brand et al. (1997),  represents predicted results from the model for nin = 1.7, nex = 5.0, S = 0.01.

5. Discussion In order to understand particle transport and deposition processes in the airways, aerosol bolus dispersion is simulated using the same breathing conditions and data analysis as the measurements by Brand et al. (1997), which have been widely used in many model simulations as a reference (Darquenne et al., 1997; Darquenne & Prisk, 2003; Goo & Kim, 2001; Kohlhaufl et al., 1999, 2000; Lee & Lee, 2002; Sarangapani & Wexler, 1999, 2000). The parameters optimized to the measurement data in this model allow us to predict bolus dispersion in terms of standard deviation, MS, skewness and deposition. It is instructive to examine the sensitivity of the model to small changes in the identified parameter values. Sensitivity to the scaling factor, S, in the cumulative mixing intensity was explored. The results for S = 0, 0.005, 0.01, and 0.015 were compared with measurement data. Both MS and standard deviation increase with increasing scaling factor whereas skewness decreases with increasing scaling factor. The higher mixing intensity increases bolus dispersion and flattens the MS dependence on the mean penetration volume. The approximate contribution of asymmetry in the particle transport profiles and chaotic advection due to secondary motion to dispersion can be estimated by comparing the results with a zero cumulative mixing intensity (S = 0) to the observation. For MS, the result with S = 0 is similar to the observations for volumetric lung depths less than 200 ml, showing that the asymmetric particle transport profiles between inhalation and exhalation contribute to the measured dispersion more than cumulative mixing intensity in the upper airways, while mixing intensity due to secondary flow on exhalation plays a more significant role to the overall measured dispersion at deeper volumetric lung depth. For standard deviation, the overall influence of cumulative mixing contributes almost 70% of the total dispersion at shallow lung depths (∼ 100 ml) and 45% at deeper lung depths. The skewness for S = 0 significantly overestimated the measurements. Therefore, it is not possible to predict bolus dispersion properly with only the mean flow pattern differences between inhalation and exhalation. We also examined how each bolus shape parameter changes in response to differing inspiratory effective particle transport exponent. All three bolus descriptors agree well with the measurements when nin is in the range of 1.7 and 2.0. The MS tends to decrease with decreasing nin , whereas standard deviation and skewness increase. This indicates that bolus mode is shifted more proximally and bolus is more dispersed as the inspiratory effective particle transport profile becomes sharper. Furthermore, the effect of changing the expiratory exponent is explored. Parameters obtained from the values between 3.0 and 8.0 for nex are compared well with the measurements. As the expiratory effective transport profile becomes blunter with decreasing nex , it freezes the developed inspiratory particle transport profile so that the bolus is more dispersed resulting in increasing standard deviation and skewness while dramatically decreasing MS as penetration volume increases.

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

243

During aerosol bolus dispersion experiments, expired particle concentration is measured at the mouth after inhaling a small portion of aerosol bolus to a given lung depth. Measured dispersion of the recovered bolus depends on the structure of the lungs, breathing maneuver, and experimental apparatus. The intersubject variability in airway morphology (airway diameters, lengths, and branching angles) and different breathing patterns may add variability to dispersion and deposition fractions even within a given airway generation. The analysis of dependencies between aerosol bolus dispersion and lung function parameters showed a positive relationship between total lung capacity and dispersion. As total lung capacity increases, bolus dispersion also increases because of more amount of air available for mixing. Since we employed Weibel’s morphometric data (Weibel, 1963; Haefeli-Bleuer & Weibel, 1988), this may also be responsible for some of the difference from the median values of experimental data. Another potential explanation for a large scatter in the measurements is that the measurements are challenging especially to maintain steady flow rate and the desired volumetric lung depth. In Fig. 8, the total deposition of ∼ 0.5 m particles is less than 20% of the total inhaled particles for the shallower breaths (∼ 383 ml of volumetric lung depth), and only 30% for the deeper breaths supporting the contention that particle deposition plays a minor role for ∼ 0.5 m particle bolus dispersion. Thus 0.5 m aerosol bolus dispersion is mostly attributed to convective mixing and not to other intrinsic particle motions (Darquenne et al., 1997; Heyder et al., 1988; Schulz et al., 1996; Taulbee et al., 1978). In summary, we developed a semi-empirical model for bolus dispersion and identified the effective inspiratory and expiratory particle transport profiles averaged over the entire tidal volume cycle, and the intensity of secondary motions as a function of the lung depth by simulating the bolus dispersion data of Brand et al. (1997). Although the simple power-law model in the present study does not attempt to predict the details of local complex air flow patterns in the respiratory airways, the characteristics of particle dispersion, which are more related to the inhaled volume than elapsed time, is well captured, so the model simulates the overall behavior of bolus transport in the lung. The predicted effective particle transport profile on inhalation is near parabolic (nin = 1.7), whereas the expiratory one is relatively blunt (nex = 5.0). As a result, a bolus is stretched during inhalation so that the bolus tip reaches a volumetric depth more than twice the inhaled volume. This shows that particles penetrate to more distal volumes in the airways than the inhaled volume would suggest. Since the exponent of the optimized expiratory transport profile is 5 instead of a value that is closer to infinity (i.e., plug flow), the transport profile is not representative of complete radial mixing over the lumen during exhalation as previous researchers had assumed (Sarangapani & Wexler, 1999, 2000). Local mixing intensities in the deep lung are very small due to low Re, but significant mixing due to secondary flows occurs mainly in the proximal airways during exhalation. Air flow in the lung may be characterized by streaming and mixing, which convectively transports particles from the inhaled air to the residual volume in the lung. On a single breath, inhaled particles rarely travel to the pulmonary regions. However, during subsequent tides retained particles are transported much deeper into the pulmonary region. In an accompanying paper of this work (Park & Wexler, 2006), we apply effective particle transport profiles and cumulative mixing intensities obtained from this work to predict pulmonary particle transport and deposition during multiple breaths and compare that work to the observations of Davies et al. (1972). Therefore, this semi-empirical single-breath model provides a framework for understanding particle transport and deposition during multiple breaths. Acknowledgments The authors wish to acknowledge helpful suggestions from Drs. Akira Tsuda, Frank Henry, and Bahman Asgharian. Research described in this article was supported by Philip Morris USA Inc. and Philip Morris International. References Anderson, P. J., Hardy, K. G., Gann, L. P., Cole, R., & Hiller, F. C. (1994). Detection of small airway dysfunction in asymptomatic smokers using aerosol bolus behavior. American Journal of Respiratory and Critical Care Medicine, 150, 995–1001. Brand, P., Rieger, C., Schulz, H., Beinert, T., & Heyder, J. (1997). Aerosol bolus dispersion in healthy subjects. European Respiratory Journal, 10, 460–467. Briant, J. K., & Lippmann, M. (1993). Aerosol bolus transport through a hollow airway cast by steady flow in different gases. Aerosol Science and Technology, 19, 27–39. Brown, J. S., Gerrity, T. R., Bennett, W. D., Kim, C. S., & House, D. E. (1995). Dispersion of aerosol boluses in the human lung—dependence on lung-volume, bolus volume, and gender. Journal of Applied Physiology, 79, 1787–1795.

244

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

Cheng, Y. S. (2003). Aerosol deposition in the extrathoracic region. Aerosol Science and Technology, 37, 659–671. Cheng, Y. S., Zhou, Y., & Chen, B. T. (1999). Particle deposition in a cast of human oral airways. Aerosol Science and Technology, 31, 286–300. Comer, J. K., Kleinstreuer, C., & Kim, C. S. (2001). Flow structures and particle deposition patterns in double-bifurcation airway models. Part 2. Aerosol transport and deposition. Journal of Fluid Mechanics, 435, 55–80. Comer, J. K., Kleinstreuer, C., & Zhang, Z. (2001). Flow structures and particle deposition patterns in double-bifurcation airway models. Part 1. Air flow fields. Journal of Fluid Mechanics, 435, 25–54. Darquenne, C. (2001). A realistic two-dimensional model of aerosol transport and deposition in the alveolar zone of the human lung. Journal of Aerosol Science, 32, 1161–1174. Darquenne, C., Brand, P., Heyder, J., & Paiva, M. (1997). Aerosol dispersion in human lung: comparison between numerical simulations and experiments for bolus tests. Journal of Applied Physiology, 83, 966–974. Darquenne, C., & Paiva, M. (1994). One-dimensional simulation of aerosol transport and deposition in the human lung. Journal of Applied Physiology, 77, 2889–2898. Darquenne, C., & Prisk, G. K. (2003). Effect of gravitational sedimentation on simulated aerosol dispersion in the human acinus. Journal of Aerosol Science, 34, 405–418. Davies, C. N., Heyder, J., & Subba Ramu, M. C. (1972). Breathing of half-micron aerosols I. Experimental. Journal of Applied Physiology, 32, 591–600. Fresconi, F. E., Wexler, A. S., & Prasad, A. K. (2003). Expiration flow in a symmetric bifurcation. Experiments in Fluids, 35, 493–501. Goo, J., & Kim, C. S. (2001). Analysis of aerosol bolus dispersion in a cyclic tube flow by finite element method. Aerosol Science and Technology, 34, 321–331. Grotberg, J. B. (2001). Respiratory fluid mechanics and transport processes. Annual Review of Biomedical Engineering, 3, 421–457. Haber, S., Yitzhak, D., & Tsuda, A. (2003). Gravitational deposition in a rhythmically expanding and contracting alveolus. Journal of Applied Physiology, 95, 657–671. Haefeli-Bleuer, D., & Weibel, E. R. (1988). Morphometry of the human pulmonary acinus. The Anatomical Record, 220, 401–414. Hammersley, J. R., & Olson, D. E. (1992). Physical models of the smaller pulmonary airways. Journal of Applied Physiology, 72, 2402–2414. Haselton, F. R., & Scherer, P. W. (1982). Flow visualization of steady streaming in oscillatory flow through a bifurcating tube. Journal of Fluid Mechanics, 123, 315–333. Henderson, Y., Chillingworth, F. P., & Whitney, J. L. (1915). The respiratory dead space. American Journal of Physiology, 38, 1–19. Henry, F. S., Butler, J. P., & Tsuda, A. (2002). Kinematically irreversible acinar flow: A departure from classical dispersive aerosol transport theories. Journal of Applied Physiology, 92, 835–845. Heyder, J., Blanchard, J. D., Feldman, H. A., & Brain, J. D. (1988). Convective mixing in human respiratory-tract—estimates with aerosol boli. Journal of Applied Physiology, 64, 1273–1278. Horsfield, K., Dart, G., Olson, D. E., Filley, G. F., & Cumming, G. (1971). Models of the human bronchial tree. Journal of Applied Physiology, 31, 207–217. Ingham, D. B. (1975). Diffusion of aerosols from a stream flowing through a cylindrical tube. Journal of Aerosol Science, 6, 125–132. Kim, C. S., Iglesias, A. J., & Garcia, L. (1989a). Deposition of inhaled particles in bifurcating airway models: I. Inspiratory deposition. Journal of Aerosol Medicine, 2, 1–14. Kim, C. S., Iglesias, A. J., & Garcia, L. (1989b). Deposition of inhaled particles in bifurcating airway models: II. Expiratory deposition. Journal of Aerosol Medicine, 2, 15–27. Kohlhaufl, M., Brand, P., Rock, C., Radons, T., Scheuch, C., Meyer, T. et al. (1999). Noninvasive diagnosis of emphysema—Aerosol morphometry and aerosol bolus dispersion in comparison to HRCT. American Journal of Respiratory and Critical Care Medicine, 160, 913–918. Kohlhaufl, M., Brand, P., Scheuch, G., Meyer, T., Schulz, H., Haussinger, K. et al. (2000). Aerosol morphometry and aerosol bolus dispersion in patients with CT-determined combined pulmonary emphysema and lung fibrosis. Journal of Aerosol Medicine, 13, 117–124. Lee, D. Y., & Lee, J. W. (2002). Dispersion of aerosol bolus during one respiration cycle in a model of lung airways. Journal of Aerosol Science, 33, 1219–1234. Park, S.S. & Wexler, A.S. (2006). Particle deposition in the pulmonary region of the human lung: Multiple breath aerosol transport and deposition. Journal of Aerosol Science, submitted for publication. Pich, J. (1972). Theories of gravitational deposition of particles from laminar flows in channel. Journal of Aerosol Science, 6, 351–361. Ramuzat, A. & Riethmuller, M.L. (2002). PIV investigations of oscillating flows within a 3D lung Multiple Bifurcation model. 11th International Symposium on Applications of Laser Techniques to Fluid Flows. Lisbon, Portugal, June 13–16, 19-1. Rosenthal, F. S., Blanchard, J. D., & Anderson, P. J. (1992). Aerosol bolus dispersion and convective mixing in human and dog lungs and physical models. Journal of Applied Physiology, 73, 862–873. Sarangapani, R., & Wexler, A. S. (1999). Modeling aerosol bolus dispersion in human airways. Journal of Aerosol Science, 30, 1345–1362. Sarangapani, R., & Wexler, A. S. (2000). The role of dispersion in particle deposition in human airways. Toxicological Sciences, 54, 229–236. Scherer, P. W., & Haselton, F. R. (1982). Convective exchange in oscillatory flow through bronchial-tree models. Journal of Applied Physiology, 53, 1023–1033. Scheuch, G., & Stahlhofen, W. (1994). Aerosol dispersion in human airways during one breathing cycle: the dependence of the aerosol penetration. Journal of Aerosol Science, 25, S553–S554. Scheuch, G., & Stahlhofen, W. (1991). Effect of heart-rate on aerosol recovery and dispersion in human conducting airways after periods of breathholding. Experimental Lung Research, 17, 763–787. Schreider, J. P., & Raabe, O. G. (1981). Structure of the human respiratory acinus. American Journal of Anatomy, 162, 221–232. Schroter, R. C., & Sudlow, M. F. (1969). Flow patterns in models of the human bronchial airways. Respiration Physiology, 7, 341–355. Schulz, H., Schulz, A., Brand, P., Tuch, T., von Mutius, E., Erdl, R. et al. (1995). Aerosol bolus dispersion and effective airway diameters in mildly asthmatic children. European Respiratory Journal, 8, 566–573.

S.S. Park, A.S. Wexler / Aerosol Science 38 (2007) 228 – 245

245

Schulz, H., Schulz, A., & Heyder, J. (1996). Influence of intrinsic particle properties on the assessment of convective gas transport by aerosol bolus technique. Experimental Lung Research, 22, 393–407. Schulz, A., Tuch, T., Brand, P., Schulz, H., Erdl, R., Vonmutius, E. et al. (1994). Aerosol bolus dispersion in the respiratory-tract of children. Experimental Lung Research, 20, 119–130. Siekmeier, R., Schiller-Scotland, C. F., Kronenberger, H., & Gebhart, J. (1992). Use of aerosol pulse parameters as a sensitive discriminator between middle-aged healthy smokers and nonsmokers. Journal of Aerosol Science, 23, S1:S483–S486. Simone, A. F., & Ultman, J. S. (1982). Longitudinal mixing by the human larynx. Respiration Physiology, 49, 187–203. Taulbee, D. B., & Yu, C. P. (1975). Theory of aerosol deposition in human respiratory-tract. Journal of Applied Physiology, 38, 77–85. Taulbee, D. B., Yu, C. P., & Heyder, J. (1978). Aerosol transport in human lung from analysis of single breaths. Journal of Applied Physiology, 44, 803–812. Ultman, J. S. (1985). Gas transport in the conducting airways. In: L. A. Engel, & M. Paiva (Eds.), Gas Mixing and Distribution in the Lung (pp. 63–136). New York: Dekker. Ultman, J. S., & Thomas, M. W. (1979). Longitudinal mixing in pulmonary airways—comparison of inspiration and expiration. Journal of Applied Physiology, 46, 799–805. Verbanck, S., Schuermans, D., & Vincken, W. (2001). Saline aerosol bolus dispersion. I. The effect of acinar airway alteration. Journal of Applied Physiology, 90, 1754–1762. Weibel, E. R. (1963). Morphometry of the human lung. New York: Academic Press. Yu, C. P. (1978). Exact analysis of aerosol deposition during steady breathing. Powder Technology, 21, 55–62. Zhang, Z., & Kleinstreuer, C. (2002). Transient airflow structures and particle transport in a sequentially branching lung airway model. Physics of Fluids, 14, 862–880. Zhao, Y., Brunskill, C. T., & Lieber, B. B. (1997). Inspiratory and expiratory steady flow analysis in a model symmetrically bifurcating airway. Journal of Biomechanical Engineering, 119, 52–58.