## Abstract

It is largely acknowledged that inhaled particles ranging from 0.001 to 10 μm are able to reach and deposit in the alveolated regions of the lungs. To date, however, the bulk of numerical studies have focused mainly on micrometer-sized particles whose transport kinematics are governed by convection and sedimentation, thereby capturing only a small fraction of the wider range of aerosols leading to acinar deposition. Too little is still known about the local acinar transport dynamics of inhaled (ultra)fine particles affected by diffusion and convection. Our study aims to fill this gap by numerically simulating the transport characteristics of particle sizes spanning three orders of magnitude (0.01-5 μm) covering diffusive, convective, and gravitational aerosol motion across a multigenerational acinar network. By characterizing the deposition patterns as a function of particle size, we find that submicrometer particles [ (0.1 μm)] reach deep into the acinar structure and are prone to deposit near alveolar openings; meanwhile, other particle sizes are restricted to accessing alveolar cavities in proximal generations. Our findings underline that a precise understanding of acinar aerosol transport, and ultrafine particles in particular, is contingent upon resolving the complex convective-diffusive interplay in determining their irreversible kinematics and local deposition sites.

- inhaled particles
- particle deposition
- particle transport
- pulmonary acinus

aerosol transport and deposition in the deep pulmonary acinar regions is paramount when considering inhalation therapy or health risks associated with occupational and environmental hazards (7, 32, 38, 48). In particular, the determination of local aerosol deposition patterns has ramifications in the context of both systemic delivery and biological homeostasis (66). Following experimental deposition studies in humans using radiolabeled particles measuring the amount of activity (i.e., clearance) remaining over time [e.g., see comprehensive review by Stahlhofen et al. (57)], it is widely acknowledged that aerosols with particle diameters in the broad range *d*_{p} = 0.001–10 μm have the ability to reach the acinar regions of the lungs (23, 24, 28). Yet, such experimental studies have been mostly limited to coarse estimations only of regional (e.g., alveolar) deposition and thus fall short of uncovering local deposition densities (i.e., surface dose) as well as temporal particle dynamics. Such knowledge is relevant whether aerosols deposit preferentially in proximal or distal regions of the acinus (51) or on the tips of alveolar septa rather than inside alveolar cavities (74), as shown from histology.

Due to submillimeter structures of the distal respiratory regions [ (10^{−4} m)], in and ex vivo experiments are still limited in both spatial and temporal resolution (58); thus classic experimental inhalation studies fall short in capturing deposition timescales and dynamic particle trajectories in the acinus. As an alternative, scaled-up experimental setups have been sought and are able to match dynamic similarity of acinar flows yet typically unable to match simultaneously similarity conditions for particle transport (15, 44). In turn, numerical simulations have proven to be an attractive strategy to shed some light on both acinar flow characteristics and acinar aerosol deposition (9, 21, 22, 35, 42, 59, 60). Most frequently, computational studies have modeled simple, geometric alveolar structures connected to an acinar duct (18, 22, 26, 36, 59, 63, 65), representing an isolated alveolus or a single acinar airway branch. In an attempt to expand such numerical models, recent efforts have succeeded in generating acinar models that integrate hundreds of alveoli spanning multiple bifurcating generations (31, 42, 43, 60). However, even these models only mimic a very small fraction of a single acinus (or 1/8 subacinus) (20). Recent advances in the resolution of imaging modalities (e.g., microcomputed tomography, etc.) have led to some of the first computational fluid dynamics (CFD) simulations using three-dimensional (3D) anatomically reconstructed acinar geometries (22, 54, 61). To date, however, these studies have been limited to terminal sacs or few alveoli surrounding an acinar duct.

Whether two-dimensional (2D) or 3D geometries of single isolated alveoli (19, 67, 64), acinar ducts (8, 21, 26, 35, 36, 63), or bifurcating acinar airways (6, 10, 43, 60), numerical studies have largely focused on the kinematics of inhaled micrometer-sized particles predominantly governed by (26, 31, 43, 60); yet, this range of sizes represents only a small fraction of the aerosols capable of reaching the acinus (23, 24). For finer submicrometer particles, their motion has long been thought to follow convective airflows (65), and thus most studies have analyzed trajectories of passive tracer particles subject to acinar flows only (5, 10, 22, 25, 53, 68, 69). However, these fail to capture intrinsic particle diffusion (i.e., Brownian motion) that is anticipated to become a critical transport mechanism in the acinus (58), in particular for ultrafine particles (UFPs) with diameters <100 nm. To the best of our knowledge, there are few if any studies that have modeled Brownian motion of inhaled particles other than 2D models limited to a terminal alveolar sac (4, 9), an acinar duct (63), or a single alveolus (39). Hence, there remains a large gap in providing a detailed quantitative mapping of acinar aerosol deposition and specifically in resolving the complex interplay between diffusive and convective mechanisms of submicrometer particles reaching the acinar region.

The present study aims to shed some light on the dynamics of particle transport and deposition in the acinar region of the lungs by resolving some of the complex dynamics that arise between convective and diffusive mechanisms for submicrometer particles and, in particular, UFPs. In an effort to provide a quantitative mapping of acinar particle dynamics, we investigate the transport characteristics of particle sizes spanning three orders of magnitude (*d*_{p} = 0.01-5 μm) across a multigenerational acinar network. By characterizing aerosol deposition using physiologically driven metrics including time for deposition, local deposition density, penetration depth, and alveolar vs. ductal deposition, our study aims to deliver the most comprehensive survey of local acinar aerosol dynamics to date.

## METHODS

#### Model geometry.

Morphometrical studies (41, 52), including more recent 3D imaging-based ones (37, 65, 71), have shown that alveoli have generally a polyhedral-like shape with single flat interalveolar septa separating adjacent alveoli; namely, this space-filling polyhedral shape is comparable to honeycombs or soap foam (72). In accordance with recent CFD studies (26, 31, 60), the seminal polyhedral alveolar duct structure was adopted to capture the irregular and complex acinar geometry (16).

Here, the air space is constructed from units of truncated octahedra, or 14-sided polyhedra, that consist of 6 squared and 8 hexagonal faces. The different polyhedra are then connected to form a space-filling air space by removing some of the interconnecting faces and thus forming alveolar openings or ducts. Since the entire acinar air space is constructed with repeating polyhedral units, the characteristic dimensions of acinar ducts and alveoli remain fixed across all acinar generations; this restriction is faithful for most distal acinar generations (26) but does not necessarily capture single, isolated alveoli populating the entrance of the acinar space past terminal bronchioles (20). Dimensions of our computational model are taken from morphometry (20), such that the sleeve diameter of an acinar duct at any generation is held constant (≈575 μm) at functional residual capacity (FRC). Similarly, ductal lengths of proximal and distal acinar generations span ∼375 to 563 μm, where the base polyhedral unit allows for a ductal length to be formed by two or three polyhedra (60).

Following the above, we constructed an asymmetric bifurcating acinar airway tree representative of distal acinar generations, as shown in Fig. 1*A*. The acinar model spans a maximum of six airway generations (schematically shown in Fig. 1*B*) and is comprised of a total of 343 polyhedra (i.e., truncated octahedra) forming 277 individual alveoli. Although this computational model represents one of the most extensive acinar models developed to date (31, 42, 43, 60), its volume represents nevertheless only a small fraction of a characteristic 1/8 subacinus; ∼7.7% of a subacinar volume of 22.8 mm^{3} (20).

#### Breathing motion.

Wall breathing expansions and contractions of the acinar domain are modeled with an idealized sinusoidal motion. While temporal asynchrony and heterogeneous wall motions as a result of a geometric hysteresis in lung expansion are known to exist (25, 29, 45), the principal mode of breathing is acknowledged to remain largely self-similar (2, 17). Mathematically, the motion of the computational domain is thus mimicked with (26, 35, 59, 60, 64, 70) (1)

where **x** is the local position (in 3D) of the wall boundary (note the coordinate system in Fig. 1), **x**_{0} is the position at *t* = 0, β is the length-scale expansion factor, and ω = 2π/τ is the (angular) breathing frequency where τ corresponds to the breathing period. Under self-similar breathing conditions, the length-scale expansion factor is given by β = (1 + Vt/FRC)^{1/3} − 1, where Vt is the tidal volume. Following this modeling approach, both fluid and particle motions (see below) arise solely as a result of the distensions of the acinar model. Here, we impose physiological breathing conditions mimicking tidal breathing of an average human adult at rest with Vt = 0.5 l, FRC = 3 l, and *T* = 3 s, such that β = 0.053 (26, 60). These breathing conditions result in a maximal surface expansion of ΔS = 11.2% relative to FRC, with a corresponding surface-to-volume ratio (S/V) of 24 mm^{−1} that remains approximately constant throughout the breathing cycle. Such observation is in general agreement with recent measurements in reconstructed acini obtained from microcomputed tomography in mice where the S/V ratio was not significantly altered during acinar expansion (37), with S/V = 47 mm^{−1} for volumes close to total lung capacity (TLC). While we find a decreased S/V ratio in our space-filling acinar domain, this difference may be attributed to the larger individual volumes and structures of alveoli in human lungs relative to mice.

#### Governing equations of fluid motion.

In the acinar region of the lungs, respiratory flows can be characterized as time dependent, incompressible, single phase, and Newtonian (58). Hence, air motion can be described using the continuum and momentum equations in an arbitrary Lagrangian-Eulerian (ALE) framework (11, 30) (2) (3)

where **u**_{mesh} is the velocity of the moving mesh given by **u**_{mesh} = *d***x**(*t*)/*dt*, following *Eq. 1*. Here, *t* is the time, **u** is the velocity field, p is the pressure, ρ_{f} is the density, and μ is the dynamic viscosity of the fluid (i.e., air). Note that for the air phase, body forces, such as gravity, are negligible (58). For unsteady acinar flows, two characteristic dimensionless numbers arise; namely the Reynolds and Womersley number. The Reynolds number represents the ratio of inertial to viscous forces and is defined as *Re* = *DU*/ν, where *D* is the equivalent ductal diameter, *U* is the characteristic velocity (defined here as the average inlet velocity at peak inhalation), and ν = μ_{f}/ρ_{f} is the fluid's kinematic viscosity. In most distal acinar generations, Reynolds numbers are typically on the order of *Re* ≈ (10^{−3}) describing laminar or creeping flow conditions (58). Here, maximal Reynolds numbers are *Re* < 0.3 at peak inhalation at the inlet of the domain. Similarly, the Womersley number, defined as α = *D*(ω/ν)^{1/2}, represents the ratio of unsteady to viscous forces (73) and is captured by subunity values under normal breathing conditions (α < 0.1). In turn, acinar flow patterns are anticipated to remain quasisteady over the inhalation/exhalation cycle for self-similar breathing conditions (49, 58, 65).

Computational simulations of the air phase were carried out using the open source C++ library OpenFOAM (Open Source Field Operation and Manipulation, Version 2.1.1). As recently described (25, 26), the Navier-Stokes equations were discretized using the finite-volume-method (FVM) in combination with a merged PISO-SIMPLE-algorithm for the pressure velocity coupling (27, 46). This fluid solver was coupled one way with a DEM-particle tracking (discrete-element-method) to model aerosol transport (see below). Spatial and temporal discretization were approximated using second-order schemes to ensure consistency and accuracy. For temporal discretization, the Crank-Nicholson method was implemented, combined with a Courant-Friedrichs-Lewy (CFL) condition for adaptive time stepping. The resulting time step sizes span 10^{−5} < Δ*t*/τ < 3 × 10^{−4}, depending on the instant in the breathing cycle.

To ensure independence of the computational grid used, a mesh convergence study was performed. Fluid velocities as well as deposition statistics for *d*_{p} = 0.5 μm (i.e., such particles are the most sensitive to convection only, see results) were analyzed and compared for four different mesh sizes ranging from 0.7 to 5.7 M tetrahedral cells. Overall, a mesh size of ∼1.5 M cells was found to approximate well acinar airflows. In line with our previous studies in isolated alveolated ducts (25, 26), the final mesh was locally refined in regions featuring high-velocity gradients, including proximal acinar generations as well as in the vicinity of the alveolar openings and at the tips of the alveolar septa. Differences in peak velocities in the most proximal acinar generation (*generation 0*, at peak inhalation, i.e., *t* = 0.25τ), where convective flows are strongest, were found to be <3% compared with the maximum mesh size investigated (i.e., 5.7 M cells). Additionally, particle deposition efficiency, i.e., number of deposited particles vs. total number of injected particles, differed by <4% after one breathing cycle.

#### Particle transport.

In the absence of hygroscopic growth and electrostatic forces, aerosol motion in the acinar region is governed by three fundamental transport mechanisms (58): convection, sedimentation, and diffusion. Restricting our attention to solid spherical particles, the linear momentum of a single aerosol is described by the differential equation (4)

where **u**_{p} is the particle velocity and *m*_{p} the particle mass. The left-hand side of *Eq. 4* expresses changes in a particle's momentum over time (i.e., acceleration), whereas the right-hand side corresponds to the sum of forces acting on the particle, namely fluid drag (**F**_{D}), gravitational sedimentation (**F**_{g}), and Brownian motion (**F**_{Br}), respectively. For a spherical particle, the drag force is defined as
(5)

where subscripts f and p correspond to fluid and particle properties, respectively. The particle diameter is *d*_{p} such that *m*_{p} = ρ_{p}π*d*_{p}^{3}/6. The parameter *C*_{C} corresponds to the Cunningham slip correction factor given by (14, 24)
(6)

where λ is the mean free path for air (λ ≈ 69 nm under ambient conditions). The net gravitational force is defined as (7)

where the gravitational acceleration is constant with **g** = (0 − 9.81 0)^{T} m/s^{2} (note the arbitrarily chosen orientation of the model in Fig. 1). Since recent studies have demonstrated that deposition efficiencies (and corresponding deposition time scales) are not affected by the direction of gravity in space-filling acinar domains (31), we limit our simulations to a single gravitational orientation. Finally, the stochastic force due to a particle's Brownian motion **F**_{Br} is expressed as (4, 40)
(8)

where **ñ** is a random vector following a normal Gaussian distribution with zero mean and unit variance. Here, the constant *S*_{0} corresponds to the amplitude of the white noise process, defined as
(9)

where *k* is the Stefan-Boltzmann constant (*k* = 1.38·10^{−23} J/K) and T is the absolute temperature in Kelvins. To first validate our stochastic model for Brownian motion, the root-mean-squared (RMS) displacement for an ensemble of 10,000 particles (*d*_{p} = 100 nm) injected into an infinite quiescent medium was compared with Einstein's theory (13). Our results indicate an error in RMS displacements of 2.7% after 10 s (see appendix a).

We implemented the given forces (*Eq. 4*) into our custom solver in OpenFOAM; note that the stochastic force (*Eqs. 8* and *9*) as well as the corrected drag force (*Eqs. 5* and *6*) were implemented in the kinematicCloud-library (47). *Equation 4* was solved using an implicit Euler scheme to ensure stability. Approximately 20,000 aerosols are injected as a bolus at the beginning of the first inhalation (*t* = 0), where particles are seeded homogeneously across a unit polyhedron at the inlet of the most proximal generation (26). For validation of convection-sedimentation transport, particle trajectories in a straight pipe were analyzed as previously described (26, 33, 60); here, the average error for the deposition efficiency over time was found to be <1% compared with the analytical solution.

To characterize spatial deposition patterns and quantify deposition “hotspots” in the acinar domain, deposition densities/concentrations were evaluated. Namely, we counted the number of neighboring particles deposited within a defined radius for any given particle, where the radius was chosen as 100 μm, approximately representative of the characteristic length scale of the face of an alveolar wall. As such, our methodology is based on counting particles, *n*, within a characteristic cell (or control volume, V_{Ch}).

#### Nondimensional particle numbers.

To cover a comprehensive range of particle sizes acknowledged to reach the acinar region (24, 28), we investigated six different particle sizes spanning *d*_{p} = 10 nm to 5 μm, representing a span over three orders of magnitude. Depending on the aerosol diameter, different transport mechanisms are known to influence, if not determine, aerosol trajectories and are captured by three dimensionless particle numbers (58), namely the Stokes, gravity, and (particle) Peclét number. Briefly, the Stokes number (*Stk*) represents the ratio of inertial to drag forces experienced by a particle and typically defined as *Stk* = ρ_{p}*d*_{p}^{2}ω*C*_{C}/18μ under oscillatory flow conditions (58). For all particle diameters considered, we note that *Stk* << 1 (see Table 1) such that aerosol deposition via inertial impaction is negligible in the acinus. Correspondingly, the Gravity number is defined as *H* = *g*_{p}*d*_{p}^{2}*C*_{C}/18μ*D*ω and may be interpreted as a parameter proportional to the ratio between the time scale for oscillatory flow (i.e., breathing period) and the characteristic gravitational sedimentation time (19, 60). In turn, we expect micrometersize particles to be strongly affected by sedimentation [*H* ∼ (1)] compared with UFPs (*H* << 1), as captured in Table 1.

Finally, the particle Peclét number defined as *Pe*_{p} = *d*_{p}*D*ω/*D*_{Diff,p} where *D*_{Diff,p} is the particle diffusion coefficient (see appendix a), compares an unsteady flow speed (*D*ω) to the particle's characteristic diffusive velocity scale (*D*_{Diff,p}/*D*). Values of *Pe*_{p} transition from a convection-dominated regime (*Pe*_{p} >> 1) near *d*_{p} = 5 μm to a diffusion-dominated one (*Pe*_{p} << 1) at *d*_{p} = 0.01 μm (see Table 1), with a transition from convection to diffusion [*Pe*_{p} ∼ (1)] in the region between *d*_{p} = 0.5–1 μm. These order of magnitude estimates suggest initially that the governing aerosol transport and ensuing deposition mechanisms in the acinus are expected to strongly vary with particle size.

## RESULTS

#### Particle kinematics.

While much work has been previously dedicated to the combined effects of convection and sedimentation on particle trajectories (26, 31, 42, 60), we first provide the reader with some qualitative intuition on characteristic convective-diffusive particle dynamics pertinent to submicrometer aerosols and, in particular, UFPs. For such particles, we recall that *Pe*_{p} << 1 (see Table 1) such that the role of Brownian motion is anticipated to become considerable.

Figure 2 exemplifies time-lapse particle positions over three breathing cycles for the case of *d*_{p} = 0.1 μm; here, color coding indicates deposited (blue) and remaining airborne (red) particles. At peak inhalation (*t* = 0.25τ), particles are transported deep into the acinar structure carried by relatively high flow velocities [*u*_{max} ≈ (0.1 m/s), *Re* ≈ 1] in the ductal regions. By the end of the first inhalation phase (t = 0.5τ), particles mainly undergo Brownian motion during instants when convective respiratory flows are nearly quiescent. This diffusive mechanism remains transient over a short window of time near flow reversal and is strikingly visible when looking at simulations of particle motion during breathing (see Supplemental Movies S1 to S4 for particles in the range *d*_{p} = 0.05-5 μm; Supplemental Material for this article is available online at the Journal website). Figure 2*E* shows that within a single breathing cycle, a significant proportion of particles spread and deposit over the entire acinar domain. As discussed further below, this phenomenon is characteristic for submicrometer particles in the acinar region (see Fig. 6).

To illustrate the different transport regimes (convective, diffusive and gravitational) as a function of particle size, we first focus on single particle trajectories; Fig. 3 presents particle tracks extracted in the first generation of the model. Whereas aerosol transport in ducts is mainly dominated by convective motion, the situation is distinctively different inside alveoli. There, convective respiratory flows are typically weak inside alveolar cavities and are known to exhibit recirculating (i.e., proximal generations) or more radial-like (i.e., distal generations) flow patterns where velocity magnitudes are up to several orders of magnitude smaller relative to ductal flows (26, 58, 59).

As a result, particle motion may be highly influenced by either diffusion (i.e., UFPs) or sedimentation (*d*_{p} ≥ 1 μm) and exemplified for four different particle sizes: Fig. 3*A* shows a highly diffusing *d*_{p} = 0.1-μm particle characterized by an intricate, jiggling trajectory stereotypical of diffusive behavior. At *d*_{p} = 0.5 μm (Fig. 3*B*), the particle follows much more closely the local flow streamlines near the alveolar opening, somewhat acting as a passive tracer during the beginning of the inhalation phase, before undergoing small diffusive hops to reach the center of the alveolus by the end of the breathing cycle; such latter trajectory highlights the intricate coupling that arises between convection and diffusion locally within an alveolus for a submicrometer particle. Meanwhile, a slightly larger particle (*d*_{p} = 1 μm) exhibits little Brownian motion and is instead pulled into the alveolar cavity during inhalation before exhibiting a small, yet nonnegligible gravitational pull once within the cavity, given the weak convective alveolar flow (Fig. 3*C*). Finally, heavier particles such as *d*_{p} = 5 μm (Fig. 3*D*) undergo fast deposition within less than an inhalation phase, under the sole influence of gravity; note the alignment of the particle trajectory with the gravity vector.

#### Particle deposition.

To evaluate the fate of inhaled micrometer and submicrometer aerosols, we first extract some of the underlying deposition statistics. Figure 4*A* presents deposition efficiencies (fraction) in the acinar domain as a function of time. For UFPs with diameters of *d*_{p} ≤ 0.05 μm, total deposition inside the domain is reached almost instantly (*t* < 0.1τ) due to their highly diffusive nature. Concurrently, heavy aerosols (*d*_{p} = 5 μm) deposit nearly just as fast under the fast action of gravity with *H* > 1 (Fig. 4*A*). In contrast, aerosols in the intermediate size range (0.1 μm < *d*_{p} ≤ 1 μm) exhibit much lower deposition efficiencies with long, asymptotic deposition times; namely, these particles are neither subject to strong diffusive nor gravitational mechanisms. For example, the deposition efficiency for *d*_{p} = 0.5 μm particles is the lowest (within the finite particle sizes investigated) and reaches only 60% after five breathing cycles, whereas 0.1- and 1-μm particles reach asymptotically <85% and <80%, respectively, over the same cumulative time. Note that particles not deposited have, for the main part, exited the domain within the first exhalation phase and are thus considered as exhaled (i.e., they do not reenter the model). As such, our results do not capture deposition phenomena occurring in more proximal airway generations (e.g., terminal bronchioles) and cannot be directly compared with semiempirical models for broad lung deposition data in the alveolar region (24, 28).

Figure 4*B* classifies deposited particles according to ductal and alveolar deposition. Here, ductal deposition includes deposition in the connecting ducts between the acinar generations as well as on the alveolar rings, where we define the thickness of the ring to be 12.5 μm (5% *d*_{Alveolar}) following in vivo bronchoscopy measurements in human adults (62). The ratio of alveolar to ductal deposition varies considerably according to particle size, where larger (heavy) and smaller (diffusive) particles in the spectrum are more likely to reach alveoli (i.e., inside alveolar cavities). In contrast, intermediate particle sizes mainly deposit near alveolar openings in the acinar ducts. Overall, we observe a significant spread in the total alveolar deposition fraction ranging from ∼14.4% (*d*_{p} = 0.5 μm) to 60.7% (*d*_{p} = 5 μm).

As Fig. 5 emphasizes (note the logarithmic scale), deposition density varies strongly according to particle diameter. Smaller and larger aerosols yield localized “hotspots,” in particular in the first acinar generation (Fig. 5, *A* and *F*). In contrast, particles with sizes in the midrange (*d*_{p} = 0.1–0.5 μm) are much more widely spread over the entire domain (Fig. 5, *C*–*E*), resulting in considerably lower deposition densities. We note that for a fixed number of initially injected particles the maximum deposition density spans an order of magnitude from 388 (*d*_{p} = 0.5 μm) to 3,910 particles (*d*_{p} = 5 μm) per control volume; these latter results are correlated with the local acinar location of density maxima. For example, while the deposition hotspot for *d*_{p} = 5 μm aerosols is located in *generation 0* (Fig. 5*D*), that of *d*_{p} = 0.5 μm is located in the third acinar generation (Fig. 5*F*).

Following the above results, Fig. 6 discriminates the spatial distribution of aerosol deposition according to their penetration depth; namely, the relative deposition (i.e., distribution) is shown according to acinar generation. As Fig. 6 shows, the deposition of small (*d*_{p} = 0.1 μm) and large (*d*_{p} = 5 μm) particles is strongly biased towards proximal acinar generations with a peak deposition in *generation 0* of more than 95 and 80%, respectively; it should be noted that for 5-μm particles local deposition sites within *generation 0* depend on the specific orientation of gravity. Nevertheless, these results emphasize that both highly diffusive and heavy particles are unlikely to reach the deeper, more distal acinar generations, i.e., diffusion and sedimentation act too quickly. Particles with diameters in the range of 0.1–1 μm deposit on the other hand in more distal acinar generations, with a peak located in *generation 2*. For such particles, a more uniform deposition distribution can be noted across the acinar generations, in line with the findings of Fig. 5. Namely, 0.1-μm aerosols exhibit the most homogeneous spread with a relative maximum of only 25.8%. Overall, however, deposition in the most distal acinar generations (i.e., *generations 4* and *5*) remains overwhelmingly low (<<10%) if not absent, across all particle sizes. This is a direct result of the finite maximal ductal path length of the convective tidal front under quiet Vt conditions (i.e., a deep inhalation maneuver is anticipated to carry the bolus front more distally).

## DISCUSSION

The majority of previous studies have been often based on the underlying assumption that inhaled submicrometer aerosols are mainly affected by convection in the acinar region, and as such (ultra)fine particles act mainly as passive flow tracers, in the absence of Brownian diffusion (10, 22, 25, 53, 54); such assumption has often followed from considerations limited to the relatively large values of the particle Peclét number (*Pe*_{p} >> 1) based on average (acinar) ductal velocities. However, estimations based on transport in ductal airways fall short of capturing localized aerosol kinematics inside alveoli or dynamics occurring at flow reversal. Few investigations have incorporated diffusive mechanisms into their particle transport model, and those that have been mostly limited to 2D domains featuring an isolated alveolar sac (4, 6) or few alveoli surrounding a duct (63). Dimensional analysis underlines that submicrometer particles, and even more so UFPs, are anticipated to be strongly affected by diffusion once they reach the acinar region (58), as highlighted by small values of *Pe*_{p} for *d*_{p} < 0.5 μm (see Table 1).

In the present study, we have attempted to shed new light on acinar particle transport covering a broad range of submicrometer- to micrometer-sized particles and thereby characterize the roles of convection, sedimentation, diffusion, and their interplay in a 3D acinar tree. In an effort to quantify the role of particle diffusion and its importance in governing submicrometer aerosol transport, simulations in the presence and absence of Brownian motion (see appendix b and Fig. B1) were performed. We find that particles as large as *d*_{p} = 0.5 μm experience a significant increase in the deposition fraction compared with a case without diffusion (Δ = 14.7% evaluated at *t* = 2τ; see Fig. B1*B*). In parallel, the influence of gravity starts to affect the deposition fraction for particles as small as *d*_{p} ≥ 0.5 μm, where we find an increase of Δ = 5.8% (at *t* = 2τ) relative to deposition in the absence of gravity (see Fig. B1*A*). Hence, the dynamics of particles larger than *d*_{p} > 0.5 μm should include gravity when considering acinar transport whereas for *d*_{p} < 1 μm, the inclusion of diffusion becomes increasingly critical to capture realistic aerosol kinematics and ensuing deposition time scales.

Unlike often considered (10, 22, 25, 53, 54), submicrometer particles (i.e., *d*_{p} = 0.5 − 1 μm) cannot be regarded as purely passive tracers; such considerations would only be valid when examining bulk transport in the ductal regions of the acinus where particle Peclét numbers are large. In particular, the role of diffusion becomes increasingly significant once particles enter individual alveoli (Fig. 3, *B* and *C*) and/or during instants near flow reversal when convection is weak or quiescent. In both circumstances, diffusion gives rise to an irreversible streamline crossing mechanism. Our simulations underline the limitations of a nondimensional analysis (Table 1) that can only provide some initial intuition on particle kinematics. Here, we see that trajectories of UFPs in alveolar cavities give rise to the significant influence of intrinsic diffusion (substituted by sedimentation in the case of heavier particles) in determining particle residence time (Fig. 4*A*) and local deposition sites (Figs. 4*B* and 6). Altogether, our numerical results exemplify the respective roles of diffusive, convective, and gravitational mechanisms on acinar aerosol transport and furthermore the complex interplay of such mechanisms (i.e., convective diffusive, convective sedimentation) leading to irreversible kinematics. While previous studies have highlighted some of the complex dynamics of micrometer-scale particles undergoing convective-sedimentation transport in acinar networks (31, 42, 43, 60), omitting Brownian motion when treating the dynamics of UFPs, and more generally submicrometer particles, leads to erroneous estimations of deposition efficiencies (see appendix b) and local alveolar deposition sites.

From a more physiological perspective, our study emphasizes the relatively low acinar deposition efficiencies for aerosols in the size range of *d*_{p} = 0.1–1 μm, with a minimum at *d*_{p} = 0.5 μm (Fig. 4*B*). This latter result is in agreement, at least qualitatively, with the deposition minimum reported for regional (alveolar) deposition from semiempirical models such as that of the International Commission on Radiological Protection (IRCP) (24, 28). Comparing our results for larger *d*_{p} = 1 μm to previous studies in alveolated airway trees, we find relatively higher deposition efficiencies compared with the work of Ma and Darquenne (43) but slightly decreased compared with early space-filling models (60), although trends remain generally similar. When considering particle residence times in the acinar model, we find that these differ significantly depending on particle size, with a maximum at *d*_{p} = 0.5 μm (Fig. 4*A*). For larger micrometer-sized particles, our study highlights a decreased residence time and increased deposition efficiency in agreement with previous works (6, 31). For UFPs with decreasing diameter, an increasing residence time was previously predicted in the absence of diffusion (31, 54). In contrast, our results support that Brownian motion substantially accelerates deposition for decreasing submicrometer particle sizes, where UFPs (*d*_{p} < 0.1 μm) with *Pe*_{p} << 1 quickly deposit within less than a single breathing cycle. Since semiempirical models such as the IRCP (24, 28) only capture deposition of UFPs across the broad alveolar region, the fast deposition dynamics observed here for *d*_{p} = 0.01 μm raise the question of whether UFPs would not already deposit in more proximal terminal bronchioles (not simulated in our domain).

When considering the spatial dispersion of inhaled particles, ductal vs. alveolar deposition is anticipated to be physiologically relevant (50, 66, 74). The deposition patterns shown in Figs. 4*B* and 6 suggest that particles that reach deeper into the acinus are likely to deposit in the duct or at the alveolar opening (i.e., alveolar ring), as suggested (9, 63). From a drug delivery perspective (12, 23, 32), our results would suggest a fundamental paradigm limiting either the depth of drug delivery in the acinus or targeting alveolar deposition (i.e, inside alveoli), bearing in mind that simulations pertain to spherical particles only. Despite such limitations, Fig. 5 is in agreement with previously reported trends (6), where we observe an increasingly heterogeneous deposition pattern for particle sizes larger than *d*_{p} > 1 μm.

Looking at in vivo measurements of inhaled UFPs, Kreyling et al. (34) found an increased translocation of *d*_{p} = 20 nm compared with *d*_{p} = 80 nm particles; following our findings (Fig. 4*B*), this could be a result of increased alveolar deposition with decreasing UFP diameter. For the entire range of particle sizes, we have shown that the penetration depth of inhaled particles (51) is strongly dependent on particle size (Fig. 6). In particular, the highest deposition efficiencies (*d*_{p} < 0.1 μm and *d*_{p} > 5 μm) were found to coincide with high local concentrations restricted to proximal acinar airways only (Fig. 5). This latter result is of interest when considering exposure to hazardous particles that presumably lead to increased tissue damage as a result of high local deposition concentrations (3, 51).

Next, we discuss some of the limitations pertinent to the scope of the present work. First, the modeled acinar geometry is based on idealized polyhedra (16, 35, 60) where our domain only represents a small fraction of a total acinus, or even subacinus; as such, our results cannot be straightforwardly extrapolated to whole lungs. Moreover, our model employs a repeating base unit to mimic the acinar structure, with fixed alveolar sizes and bifurcation angles. In reality, however, the morphology of the acinus is known to be significantly more heterogeneous (71), where acinar ducts not only become more densely populated with acinar depth but alveoli also increase slightly in size with increasing generation (20). Thus we would anticipate a minor shift towards more distal deposition for a more anatomically realistic domain with larger alveolar diameters with increasing acinar generation.

Computationally, our implemented sinusoidal and self-similar breathing motion is idealized, despite the principal self-similar mode of breathing (2, 17), and thus does not capture effects such as asymmetric expansion (25, 37, 55), heterogeneous alveolar recruitment (56), or geometric hysteresis (45); these factors are anticipated to influence significantly local deposition outcomes and enhance the irreversible nature of inhaled particle kinematics (66). Finally, we note that our aerosol bolus inhalation relies on injecting particles in synchrony with the onset of the first inhalation phase (54). While this injection is not physiologically realistic, i.e., particles are anticipated to enter the deep acinar region later in the inhalation phase after crossing through the dead space and are likely to display an inhomogeneous spatial distribution across the duct upon entering the acinus (also, particles may remain in the acinus from previous inhalations), our bolus injection method is intended to provide a comprehensive and quantitative reading of the anticipated characteristics of particle transport in such a network. We hypothesize that a more physiologically realistic injection would entail increased proximal deposition as well as decreased deposition efficiencies, due to an increased particle concentration in the center of the duct at the detriment of the duct periphery (near the walls).

By investigating a broad range of particles spanning submicrometer- to micrometer-sized aerosols, our efforts help deliver a more comprehensive understanding of the intricate coupling arising between intrinsic particle transport mechanisms and the underlying acinar airflows that ultimately governs acinar deposition outcomes. On the one hand, aerosol sizes near *d*_{p} ≈ 0.5 μm reach deepest into the acinar structure with low local deposition densities leading to rather homogeneous particle distribution biased towards ductal (or alveolar ring) deposition. These particles are retained longer in the acinus with decreased deposition efficiencies. On the other hand, particles with *d*_{p} < 0.1 μm or *d*_{p} > 1 μm show contrasting deposition characteristics with maximal deposition efficiencies that are bound to proximal generations. Overall, these results could constitute a first step towards designing inhalation therapy strategies for targeting locally the acinar region using submicrometer aerosols.

## GRANTS

This work was supported in part by the European Commission (FP7 Program) through Career Integration Grant PCIG09-GA-2011-293604 and the Israel Science Foundation Grant No. 990/12.

## DISCLOSURES

No conflicts of interest, financial or otherwise, are declared by the author(s).

## AUTHOR CONTRIBUTIONS

Author contributions: P.H. performed experiments; P.H. analyzed data; P.H. and J.S. interpreted results of experiments; P.H. prepared figures; P.H. and J.S. drafted manuscript; P.H. and J.S. edited and revised manuscript; P.H. and J.S. approved final version of manuscript; J.S. conception and design of research.

## ACKNOWLEDGMENTS

We thank R. Fishler (Technion) and Dr. J. Oakes (University of California Berkeley) for fruitful discussions.

## Appendix A

To validate our stochastic model for particle diffusion, RMS displacements of particles originating from a single position were assessed and compared with analytical predictions for particle ensembles. Theoretically, the RMS displacement of an ensemble of spherical particles under Brownian motion is given by (13) (10)

in a 3D, infinite and quiescent medium. Here, the diffusion coefficient is defined as (1) (11)

To ensure that simulations were independent of numerical constraints and statistical uncertainties, different mesh densities and time steps were examined. The particle size was held constant at *d*_{p} = 0.1 μm for a simulation time of *T* = 10 s. The maximum deviation after the total simulation time was ΔRMS = 2.7% (Fig. A1), confirming good agreement between the numerical scheme and the theoretical prediction.

## Appendix B

To investigate how the omission of gravitational or diffusional mechanisms influences acinar deposition, simulations neglecting such forces were performed in the acinar domain as a function of particle size. In particular, we determine the error (Δ) in the deposition fraction over two complete breathing cycles (*t* = 2τ) associated with neglecting each force independently (i.e., sedimentation or diffusion). Simulations were performed where we anticipate the transition to lie for particle sizes where transport regimes are captured by *Pe*_{p} ≈ 1 or *H* ≈ 1 (Table 1). In Fig. B1, we present results for the deposition fraction in the absence (and presence) of gravity (Fig. B1*A*) and diffusion (Fig. B1*B*), respectively.

In assessing the role of gravitational sedimentation (Fig. B1*A*), we find that the omission of gravity for fine *d*_{p} = 1 μm aerosols leads to important discrepancies (Δ = 36.4%). As anticipated, no significant effect is found for ultrafine particles of size *d*_{p} = 0.1 μm (Δ = 0.3%), we observe however, a small, yet nonnegligible effect beginning near *d*_{p} = 0.5 μm (Δ = 5.8%). Such results underline how the influence of gravity is tangible within our acinar domain for particle diameters in the range *d*_{p} ≥ 0.5 μm. In the case of diffusion (Fig. B1*B*), we find that the role of Brownian motion on deposition fraction is rather small for fine *d*_{p} = 1 μm particles (Δ = 3.1%). In contrast, the influence of diffusion for sub-micrometer particles such as *d*_{p} = 0.5 μm increases significantly (Δ = 14.7%), unlike previously assumed (6). This observation is in agreement with dimensionless numbers presented in Table 1, where a transition occurs in the range spanning *d*_{p} = 0.5 to 1 μm. Overall, our findings accentuate the necessity of a particle model that integrates intrinsic particle transport mechanisms (i.e., sedimentation and diffusion) to capture the anticipated kinematics of submicrometer particle in the pulmonary acinus.

- Copyright © 2015 the American Physiological Society

## REFERENCES

- 1.↵
- 2.↵
- 3.↵
- 4.↵
- 5.↵
- 6.↵
- 7.↵
- 8.↵
- 9.↵
- 10.↵
- 11.↵
- 12.↵
- 13.↵
- 14.↵
- 15.↵
- 16.↵
- 17.↵
- 18.↵
- 19.↵
- 20.↵
- 21.↵
- 22.↵
- 23.↵
- 24.↵
- 25.↵
- 26.↵
- 27.↵
- 28.↵
- 29.↵
- 30.↵
- 31.↵
- 32.↵
- 33.↵
- 34.↵
- 35.↵
- 36.↵
- 37.↵
- 38.↵
- 39.↵
- 40.↵
- 41.↵
- 42.↵
- 43.↵
- 44.↵
- 45.↵
- 46.↵
- 47.↵
- 48.↵
- 49.↵
- 50.↵
- 51.↵
- 52.↵
- 53.↵
- 54.↵
- 55.↵
- 56.↵
- 57.↵
- 58.↵
- 59.↵
- 60.↵
- 61.↵
- 62.↵
- 63.↵
- 64.↵
- 65.↵
- 66.↵
- 67.↵
- 68.↵
- 69.↵
- 70.↵
- 71.↵
- 72.↵
- 73.↵
- 74.↵