scieee AI-readable full text Open interactive document viewer

Preprint for publication: Transport and Deposition of Inhaled Fibres in a Realistic Female Airway Model: A Combined Experimental and Numerical Study

Prinz, František; Kánská, Jana; Elcner, Jakub; Hájek, Ondřej; Kummerländer, Adrian; Krause, Mathias J.; Jicha, Miroslav; Lizal, Frantisek

Full text

PREPRINT Highlights Transport and Deposition of Inhaled Fibers in a Realistic Female Airway Model: A Combined Experimental and Numerical Study Frantiˇsek Prinz, Jana K´ansk´a, Jakub Elcner, Ondˇrej H´ajek, Adrian Kummerl¨ander, Mathias J. Krause, Miroslav J´ıcha, Frantiˇsek L´ızal •ELER coupled with LBM was successfully applied in a realistic female airway model •Numerical results slightly overestimated deposition, mainly in bifurcations •Time-dependent deposition analysis provides insights for optimizing aerosol delivery •Limitations of techniques for modeling fibers in complex flows were identified •Orientation-dependent calculation is crucial for accurate deposition prediction PREPRINT Transport and Deposition of Inhaled Fibers in a Realistic Female Airway Model: A Combined Experimental and Numerical Study Frantiˇsek Prinza,∗, Jana K´ansk´aa, Jakub Elcnera, Ondˇrej H´ajeka, Adrian Kummerl¨anderb, Mathias J. Krauseb, Miroslav J´ıchaa, Frantiˇsek L´ızala aBrno University of Technology, Technicka 2896, Brno, 616 69, Czech Republic bKarlsruhe Institute of Technology, Kaiserstraße 12, Karlsruhe, 76131, Germany Abstract This study presents a combined experimental and numerical investigation of fiber transport and deposition in a realistic model of the female respiratory tract, extending to the seventh generation of branching. Numerical simulations were performed using the Euler-Lagrange Euler-Rotation (ELER) method, an efficient alternative to conventional Finite Volume Methods that benefits from explicit formulation and vast scalability, enabling fast parallelization on high-performance clusters. The ELER method was coupled with the Lattice Boltzmann Method (LBM) to simulate fiber dynamics under a realistic inspiratory flow profile. Experimental validation was conducted using an identical physical airway replica. The results demonstrated good agreement between simulations and experiments in the upper airways and trachea, with some discrepancies in the bifurcations, likely owing to the challenges of modeling complex turbulent flow with ELER. This method is more accurate than corresponding effective diameter simulations. Deposition patterns were analyzed as a function of fiber dimensions, revealing higher accuracy of the ELER method for smaller particles and confirming the tendency of higher aspect ratio fibers to penetrate deeper into the lungs. The orientation-dependent deposition mechanism was deployed, underscoring the importance of solving the actual orientations of the fibers. While advancing our understanding of fiber transport in female airways, the findings also reveal limitations in current numerical techniques, particularly in bifurcations. This study emphasizes the distinct behavior of fibrous versus spherical particles, with fibers exhibiting a greater propensity to reach deeper lung regions, which has significant implications for inhalation toxicology and drug delivery. Keywords: fiber transport, deposition, in vitro, in silico, female airway geometry, Lattice Boltzmann Method, Euler-Lagrange Euler-Rotation 1. Introduction The study of particle transport and deposition within the human respiratory system is essential for understanding various physiological and pathological processes as well as for developing effective drug delivery strategies. While considerable research has focused on male airway models, there is a growing recognition of the need to investigate female-specific airway geometries due to significant ∗Corresponding author Email address: [email protected] (Frantiˇsek Prinz) Preprint submitted to Computers in Biology and Medicine May 13, 2025 PREPRINT anatomical and physiological differences between sexes [1, 2, 3]. Females exhibit smaller airway diameters, particularly in the pharynx and velopharynx [4], differences in lung volume and breathing patterns [5], and hormonal influences that can affect airway responsiveness [6]. These differences can influence airflow patterns, particle deposition, and ultimately, the efficacy and safety of inhaled therapies [7]. To address this gap in our understanding of female-specific aerosol deposition, our study utilizes a realistic female airway model supplemented with tidal volume and a breathing pattern typical for females. Fibrous particle transport has significant implications in both safety and engineering applications. Various man-made fibers are crucial components of composite materials and filtration systems used across many industries, such as carbon fibers in aerospace composites and glass fibers in air filtration systems. Although exposure to toxic fibers such as asbestos is a severe health hazard, biodegradable fibers offer unique possibilities for targeted drug delivery [8]. Numerical simulation of the transport of non-spherical particles is more complex than of spherical particles due to the nonsymmetric shape, and hence experimental studies in vitro are indispensable. However, only a few experimental studies focusing on the deposition of fibers in human airways have been conducted. Marijnissen et al. [9] and Myojo and Takaya [10] used simplified models of a bifurcation for deposition experiments. The latter experimentally examined the deposition of fibers in Weibel’s [11] idealized geometry up to the fourth generation of airway branching, even with cyclic breathing. Su and Cheng [12] conducted experimental investigations on the deposition of carbon fibers with lengths ranging from 10 to 150 µm, using a realistic airway model extending beyond the fourth generation. At a flow rate of 15 l min−1, the majority of fibers were deposited in the upper respiratory tract, comprising the oral cavity, pharynx, and larynx, whereas 63 % passed through the replica without deposition. At higher flow rates, 43.5 l min−1and 60 l min−1, the deposition rate increased, with the fibers primarily deposited in the oropharynx and larynx. The deposition efficiency in the oral cavity was significantly lower than that in the nasal cavity. These findings agree with the results of a later simulation study by Feng et al. [13]. Subsequently, Su and Cheng [14] demonstrated that inertial impaction was the dominant deposition mechanism in the oral cavity. Furthermore, Zhou et al. [15] performed experiments using two lung models: one with the oral cavity half open and the other fully open, where the degree of opening influenced the angle of attachment to the oropharynx. Su and Cheng [16] conducted experiments using three types of materials, namely glass, TiO2, and carbon fibers, at three distinct inhalation flow rates. In contrast to their previous studies focused on carbon particles, the fibers exhibited minimal deposition within the model, irrespective of their length or inhalation flow rate. The authors attributed these results to the relatively low momentum of the lightweight fibers, which enabled them to more effectively follow airflow dynamics in the airways. They indicated that these fibers do not demonstrate a significant deposition within the airways up to the third generation of branching, thereby facilitating their penetration deeper into the lungs. In contrast, fibers with greater momentum were observed to deposit in the oral cavity, allowing for potential removal through swallowing. Bˇelka et al. [17] experimentally investigated the deposition of glass fibers in a respiratory tract replica extending from the oral cavity to the 7th bifurcation. The fiber deposition within the model was minimal, with deposition fractions of 0.7 % at a steady flow rate of 15 l min−1, 1.9 % at 30 l min−1, and 4.6% at 60 l min−1. Higher deposition occurred in the oral cavity and more complex segments. Lizal et al. [18] conducted one of the rare experimental studies focusing on the visualization of 2 PREPRINT fiber movement in the simplified model of the respiratory tract. They produced a high-speed camera recording of fibers at specific recording frames in a straight tube and upstream and downstream of the bifurcation for several flow rates. A robust statistical evaluation showed two dominant orientations, parallel and perpendicular to the streamlines, and a significantly higher probability of flips downstream of the bifurcation than those in the straight tube. In addition to experimental studies, several researchers have studied fiber transport and deposition employing computational fluid and particle dynamics (CFPD) simulations in silico. The available simulation approaches were summarized by Feng and Kleinstreuer [19]. A popular choice is the one-way coupling simplified Euler - Lagrangian approach. Here, the fibers are treated as dilute, and rotational movement is neglected. The primary influence on the trajectory comes from the computation of the shape factor CDwhere either empirical correlations (e.g., Lasso and Weidman [20], Haider and Levenspiel (H-L) [21], Tran-Cong (T-C) [22]), or theoretical formulations (Stober [23]) are used. Inthavong [24] compared the H-L and T-C models in simulating the nasal cavity at a rather low flow rate of 7.5 l min−1and showed good agreement for particles smaller than 100 µm. This study was extended by Farkas et al. [25] to a realistic male model using the same airway replica as that used by Belka et al. [17] for three stationary flow rates. Chen et al. [26] investigated a single idealized bifurcation under steady laminar flow conditions. They showed that gravity affects the deposition characteristics. Orientations of deposited fibers were obtained from micrographs. The prevalent orientation was parallel with deviations mainly in the bifurcations emphasizing the effect of interception. Supporting numerical simulations used equivalent sphere approaches. With the increase in computational power in recent decades, methods that account for rotational motion have been used accordingly. Tian et al. [27] presented the Euler-Lagrange Euler-Rotation (ELER) approach and validated it by means of the deposition fraction in a straight tube for low Reynolds numbers. In a follow-up study, this team [28] numerically investigated carbon fibers in the range of aspect ratio 1 < β < 80 on a more extensive idealized model, from the trachea to the third bifurcation, and showed different rotational behaviors when passing the bifurcations. Shanley et al. [29] deployed ELER in a realistic nasal airway model under steady laminar conditions and proposed two empirical models for pressure drop and deposition efficiency. Shachar-Berman et al. [8] investigated the total breathing maneuver mimicking the Dry Powder Inhaler (DPI) profile with a peak flow rate of 90 l min−1for an average adult in a semi-realistic model up to the 9th generation of branching, and ideal fiber dimensions for medical treatment were suggested. Li et al. made significant contributions to the ELER modeling of non-spherical particles in the respiratory tract. In their initial study [30], they simulated fiber transport in a realistic nasal cavity, demonstrating complex translational and rotational behavior. While the fibers were generally aligned with streamlines, occasional quick flips were observed. Subtle deviations in trajectory and rotation can significantly impact deposition patterns, highlighting the importance of accurate interception modeling. In a follow-up study, Li et al. [31] investigated the effects of shearinduced lift forces on non-spherical particles in a circular duct. This force caused a lateral drift, potentially dominating the deposition forces and emphasizing the importance of its consideration. In their most recent study, Li et al. [32] numerically investigated fiber transport and deposition in an extended human airway model up to the 15th generation, which is a rare scope for such studies. The evaluation of the deposition for three different aspect ratios confirmed the significance of both the aspect ratio and aerodynamic diameter in influencing the deposition curves. Kiasadegh et al. [33] compared steady and cyclic regimes in a realistic model of a 24-year3 PREPRINT old female airways down to the trachea. They emphasized the importance of transient simulation, showing significant differences in deposition fractions and penetration depth compared with steadystate simulations, especially for fibers heading deeper into the lungs. In that study, a sinusoidal inspiration profile was prescribed. Tavakol et al. [34] investigated fiber deposition in nasal cavities using the ELER method coupled with Reynolds Averaged Navier-Stokes simulations, incorporating random walk models to account for turbulent fluctuations under steady flow conditions. In a consecutive study, Abolhassantash et al. [35] deployed force formulations for non-creeping flow conditions, derived by Zastawny et al. [36] and Ouchene et al. [37], and compared them with conventional creeping flow formulations within the ELER framework for simulations in a female nasal passage under steady conditions. Greater differences between the non-creeping and creeping formulations were observed for flow rates above 20 l min−1and with increasing fiber aspect ratio. Furthermore, a comprehensive comparison of different random walk models used in conjunction with RANS simulations was conducted by Mofakham et al. [38]. Finally, recent progress in the computational modeling of fiber transport was summarized by Tian and Ahmadi [39]. Most of the aforementioned studies have been either experimental or numerical. In contrast, our combined experimental and computational research not only allows for a comparison of the total deposition in each segment of the airways but also enables a broader analysis depending on the fiber dimensions. The utilization of female realistic geometry with a realistic inspiration cycle enables focusing on the differences compared to male geometry, which has been commonly studied in the past. Building upon the previous studies of Henn et al. [40] and Prinz et al. [41], this study extends the application of the Lattice Boltzmann Method (LBM) beyond spherical particles. Specifically, the ELER approach is employed as the particle tracking algorithm and coupled with the LBM to simulate the transport and deposition of fibrous particles. This novel combination allows for a more accurate representation of fiber dynamics compared with simplified models, particularly in complex flow regimes found in bifurcations. Specifically, we aim to answer the following research questions: What is the distribution of deposited inhaled fibers within various segments of the female respiratory tract? What is the influence of fiber dimensions (length and diameter) on the deposition patterns? How do the experimental results compare to the numerical simulations using the LBMELER approach, and what are the potential sources of discrepancies? How does particle release time during the inspiration cycle affect the deposition fraction and location? The remainder of this paper is organized as follows. Section 2 describes the experimental setup and numerical methods employed, including the LBM-ELER approach and geometry used. Section 3 presents the results of both the experimental and numerical investigations, followed by a discussion comparing the two approaches and analyzing the influence of fiber dimensions and release time. Section 4 discusses the limitations of the current study. Section 5 summarizes the key findings and concludes the study. 2. Methods 2.1. Experimental setup This section describes the experimental setup used to investigate the deposition of glass fibers in a realistic female airway model. The scheme of the experimental setup is shown in Figure 1 and described in detail below. Polydisperse glass fibers with diameters ranging from 1 to 10 µm and lengths ranging from 5 to 100 µm were used in the experiments. The samples were obtained by crushing glass wool (Supafil 4 PREPRINT Output flow S13 to S22 Air intake and exhaust To breathing simulator Breathing simulator LLL LUL RLL RML RUL 10 Output filters S18 S15 S14 S13 S17 S19 S20 S16 S21 S22 Hopper Fluidized-bed type disperser Charge equilibrator Rotameter Air Glass beads overflow Figure 1: Schematic of the experimental setup used to investigate fiber deposition in a realistic female airway model. Loft, Knauf Insulation GmbH) in a crusher. The fibers were mixed with glass beads (Ballotini Impact Beads, Potter Industries Inc.) and sieved to improve the ensuing deagglomeration and dispersion of fibers within the airway replica. The mixture consisted of 2% fibers and 98% beads. The prepared mixture of glass fibers and beads was introduced into the experimental setup using a hopper, as shown in Figure 1. From the hopper, the mixture flowed into a rotary feeder, which delivered it to a fluidized-bed disperser. Within the disperser, the fibers and beads were separated, and the aerosolized fibers were then transported to a charge equilibrator (NEKR-10, Eckert and Ziegler CESIO) to neutralize their electrostatic charges. The neutralized fibers were subsequently conveyed to the airway replica through a tube connected to the air intake and exhaust system of the laboratory. This connection ensured a controlled environment and prevented the leakage of fibers into the laboratory. A realistic replica of the female airway, used here for the first time, was employed in this study. It was derived by applying a uniform linear scaling factor of 0.88 to a previously validated male airway model [42, 17]. This factor was determined through a comparative analysis of key anatomical dimensions reported in the literature for representative adult male and female respiratory tracts [43, 44, 45]. Specifically, metrics such as tracheal length and mean luminal cross-sectional area served as primary parameters considered in deriving this average scaling factor to represent overall size differences between the sexes. The replica spans from the oral cavity to the seventh generation of the lung bifurcation, and the remainder of the lung is replaced by output filters. The model was capable of simulating breathing patterns using an original breathing simulator, which reproduced the specific patterns of the five lung lobes: the right upper lobe (RUL), right middle lobe (RML), right lower lobe (RLL), left upper lobe (LUL), and left lower lobe (LLL). A model of the geometry with its segmentation is depicted in Figure 2. To enhance particle adhesion and better simulate the trapping effect of the mucus layer found in vivo, the inner surfaces of the airway replica were coated with a thin layer of silicone oil prior to each experimental run. This surface treatment is a common practice in similar in vitro deposition studies [16, 17]. 5 PREPRINT 0 1 2 3 45 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 RUL RML RLL LUL LLL G0-1 G1 G1-2 G2-3 G3-7 outlets G3-4 Figure 2: Segmentation of the realistic female airway geometry with corresponding segment numbers. Each color represents a different airway generation, as indicated in the legend. To the best of our knowledge, experimental data on fiber deposition during realistic female inspiration are not available in the literature; therefore, we conducted experiments to fill this gap. The realistic breathing profiles used in this study (see Figure 3) were derived from the data presented in the Annals of the ICRP [46]. Specifically, we scaled down the values reported in [46] by a factor of 0.8, based on a comparison of the tidal volumes for the male and female geometries in [46] and the realistic male airway model from which the current female geometry was derived [17]. The distribution of flow rates among the lung lobes was adapted from Jahani et al. [47], who analyzed four-dimensional computed tomography data of six breathing cycles averaged over one minute from healthy volunteers (50% women). The airway replica was then exposed to aerosolized fibers. After exposure, two types of samples were collected: a) samples from the output filters, which captured the fibers that passed through the entire replica, and b) samples obtained from each individual segment of the replica. To collect the fibers deposited within the replica, it was first disassembled. Each segment was then placed in a beaker and immersed in isopropanol. The beaker was placed in an ultrasonic bath to dislodge any fibers adhering to the segment walls. This process created a suspension of fibers in isopropanol. To analyze the fibers, the suspension was filtered through nitrocellulose membrane filters using a vacuum filtration pump. The filters were then dried and made transparent by placing them on a glass slide using an acetone vaporizer (QuickFix, EMS, USA). Samples from the output filters were made transparent in the same manner. Each sample was then manually evaluated under a phase-contrast microscope (Nikon Eclipse E200, Nikon, Tokyo, Japan) using a 40x objective and 6 PREPRINT 0 0.5 1 1.5 0 2 4 6 8 Time (s) Flow rate (l min−1) RUL RML RLL LUL LLL Figure 3: Flow rate profiles for each lung lobe during the inspiration cycle [46, 47]. Walton-Beckett graticule (cf. representative sample Figure A.22 in Appendix). Fiber counting was performed according to a method based on the WHO guidelines [48]. To ensure a representative count, 20 randomly selected areas were analyzed for each sample. Because of the polydispersity of the fibers, the diameter and length of each fiber were measured and recorded for further analysis. 2.2. Numerical simulation setup This section outlines the numerical methods used to simulate the transport and deposition of fibers within the female airway model. Particle transport and deposition were numerically simulated using Computational Fluid and Particle Dynamics (CFPD). The solution was divided into two parts – the fluid flow solution and particle motion solution – coupled with the velocity field data, as described below. The simulation is considered as one-way coupling, and the discrete phase represented by the particles is dilute with a very low volume fraction (approximately 10−8) and hence does not affect the flow field. 2.2.1. Fluid phase Fluid flow modeling is based on the Lattice Boltzmann Method (LBM), a mesoscopic numerical approach to transport problems based on a discretization of the Boltzmann equation [49]. Due to its algorithmic structure, the LBM is uniquely suited for highly parallel execution on state-of-theart high-performance computers [50]. Specifically, it is an efficient alternative [51] to conventional finite-volume methods, providing 32-fold performance improvements compared to OpenFOAM when fixing the numerical error in a fair comparison of an industrial reference case. In LBM, the spatial simulation domain is discretized by a regular lattice on which populations f(x,  ξ, t), describing the state of the system1, propagate along discrete velocities cifollowing a so 1f(x,  ξ, t) explicitly describes the probability of the total mass of particles in position x with microscopic velocity  ξin time t. 7 PREPRINT called collision step that relaxes the per-cell populations towards their macroscopic equilibrium distribution: fi(x +ci∆t, t + ∆t) = 1−∆t τfi(x, t) + ∆t τf(eq) i(x, t),∀i∈ ⟨0,18⟩.(1) Here, cidenotes the discrete velocity stemming from the velocity set D3Q19, f(eq) iis the equilibrium distribution function and τthe relaxation time in the Bhatnagar-Gross-Krook (BGK) collision operator resolving the particle interactions. This equation can be divided into two steps, collision and streaming. fcoll i(x, t):= 1−∆t τfi(x, t) + ∆t τf(eq) i(x, t),(2) fstr i(x +ci∆t, t + ∆t):=fcoll i(x, t).(3) Macroscopic quantities are computed using discrete velocity moments ρ(x, t) = X i fi, u(x, t) = 1 ρX i fi ξi, p(x, t) = ρ c2 s ,(4) where ρdenotes the density of the fluid, u denotes the macroscopic velocity, pdenotes the pressure, and csdenotes the speed of sound. The Chapman-Enskog expansion [52] can be used to demonstrate the convergence of this method to solutions of the incompressible Navier-Stokes equations, justifying its use in the present application. To model turbulent phenomena, a Large Eddy Simulation (LES) with the Smagorinsky subgridscale (SGS) model [53] was enabled (further implementation details are provided in [41]). The Smagorinsky constant was assigned the common value of 0.1. LES provides the time-varying largescale turbulent velocity fields, which can capture the majority (e.g., >80%) of the turbulent kinetic energy in well-resolved simulations [54]. The influence of smaller subgrid eddies on particle dispersion is commonly accounted for in Reynolds-Averaged Navier-Stokes (RANS) simulations through models such as random walk approaches [55]. For LES, while SGS turbulent dispersion effects can be modeled, their impact on the deposition of micrometer-sized particles has been reported as relatively small in certain contexts (e.g., [56, 57]. Consequently, explicit modeling of SGS particle turbulent dispersion has often not been included in LES studies of airflow and deposition in the respiratory tract ([58]). Supporting this approach, Koullapis et al. [59] demonstrated in a benchmark comparison that LES models without explicit SGS turbulent dispersion for particles achieved reasonable accuracy against experimental data for regional deposition in human airways. Based on these considerations from the literature and the common practice in the field, explicit modeling of SGS turbulent dispersion effects on particle trajectories was not included in the present study, though it remains an important area for future research. A pressure boundary condition, as proposed by Skordos [60], was imposed at the mouth inlet. At the outlets, uniform velocity boundary conditions, as proposed by Skordos [60], were prescribed, with time-dependent velocities calculated from a realistic inspiration profile (Figure 3) and the corresponding tidal volume distributions. Velocities were implemented using a linear interpolation scheme. Two consecutive inhalations were simulated with a reduced hold-up time of 0.25 s between them. Each lung lobe in our model contains two funnels, which are connected to the respective lobe pistons of the breathing simulator. The flow rate distribution between the funnels in the 8 PREPRINT Collision step Streaming step Boundary treatment New time step Moment update Flow field data postprocessing Flow initialisation Data storage Compute particle torques Compute particle forces Deposition treatment Particle data postprocessing Fiber initialisation Update angular velocity LBM ELER velocity transfer no data transfer Repeat until termination criteria are fulfilled Update orientation Update position Figure 6: Schematic of the numerical setup. After the initialization step, each iteration involves one LBM cycle and one ELER cycle. After computing the flow velocities (u) in the LBM cycle, these values are transferred to the ELER algorithm. Due to the one-way coupling, no data are transferred back from the ELER algorithm to the LBM solver. In principle, different numbers of LBM or ELER cycles can be performed within each iteration (e.g., inner iterations), as indicated by the dotted arrows. The simulation terminates when the specified criteria (e.g., simulation time or number of iterations) are met, and the data are then post-processed. proportional to the actual flow rate at that time (see Figure 3), reflecting the fact that a higher flow rate carries more particles into the airway. The numerical simulations were performed using the open-source C++ library OpenLB [70, 71, 72]. This software enables flexible and performant simulations using the LBM, benefiting from efficient parallelization and scalability of both CPUs and GPUs on high-performance computers. The suitability of OpenLB for this type of application was demonstrated in a previous paper by the authors [41], which resulted in high performance and accuracy. To simulate fiber transport and deposition, OpenLB was extended using an integrated in-house code that implemented the ELER model. The simulations were performed on the Karolina supercomputer petascale system at the IT4Innovation National Computing Center in Czechia using 1024 CPU cores distributed across multiple nodes, leveraging the high memory capacity and interconnect bandwidth of the HPC system. For comparison with experimental results, the following statistical quantities were evaluated: Deposition fraction DF = Ns Ntot ,(33) where Nsis the number of particles deposited in a given segment, and Ntot is the total number of particles entering the geometry. 15 PREPRINT The deposition efficiency DE describes the efficiency of a segment in capturing the particles. DE = Ns Ns,in ,(34) where Ns,in is the number of particles entering the segment. Volume equivalent diameter deq =3 r6V π,(35) is defined as the diameter of the spherical particle having the same volume as the fiber. Aerodynamic diameter dae =deqrρp ρ0χR (36) is defined as the diameter of a sphere with unit density (ρ0) that settles with the same terminal velocity as the fiber. χRdenotes the dynamic shape factor for a random orientation. The Stokes number, given by Stk = ρ0d2 aeu0 18µD0 ,(37) is an important dimensionless parameter for evaluating the rate of the inertial impaction. It is particularly useful for comparing data across different flow rates and geometries because it accounts for the combined effects of particle inertia, fluid viscosity, and characteristic length scales. In this equation, u0represents the characteristic velocity of the flow, and D0represents the characteristic dimension of the geometry. 2.3. Numerical verification Verification of the ELER method implementation within the LBM framework was performed using a benchmark case of laminar airflow through a circular tube. This benchmark, often with small variations, has been used for verification in several studies [27, 13, 73, 74, 75]. In this paper, the specific setup first described by Tian et al. [27] was simulated, and data from Feng et al. [13] were also used for comparison. This involved simulating laminar flow through a horizontal pipe with a diameter of 4.2 mm and a corresponding Reynolds number of Re = 169. A spheroidal particle with a minor axis of 0.5 µm and an aspect ratio of 14 was injected at a position 0.45 mm above the bottom edge of the pipe, with an initial perpendicular orientation. The particle trajectory was then tracked for 0.2 s (see Figure 7). To compare the results of our simulation with those of Tian et al. [27] and Feng et al. [13], the directional cosines between the fiber symmetry axis and the coordinate axes were evaluated. The time evolution of these quantities over 0.2 s of the simulation is shown in Figures 8(a) and 8(b). As can be seen, the fiber quickly tilts from its initial perpendicular orientation to a position parallel to the flow. After each period of 0.055 s, it undergoes a 180-degree flip. The sedimentation velocity, shown in Figure 8(c), was also compared. Gravity acts on the fiber, causing it to move in the negative y direction. The sedimentation velocity increases from zero to a terminal velocity, at which point the gravitational force is balanced by drag force. When the fiber 16 PREPRINT 0.45 mm x y z g u Figure 7: Initial condition of the fiber in the numerical verification case, showing its position and orientation relative to the airflow in the horizontal pipe. 0 0.05 0.1 0.15 0.2 −1 −0.5 0 0.5 1 Time (s) cos x(-) Tian Feng ELER (a) 0 0.05 0.1 0.15 0.2 −1 −0.5 0 0.5 1 Time (s) cos y(-) Tian Feng ELER (b) 0 0.05 0.1 0.15 0.2 −3 −2 −1 0·10−4 Time (s) vy(m s−1) Tian Feng ELER (c) Figure 8: Comparison of the time evolution of (a) the directional cosine between the symmetry axis of a fiber and the x-axis (cos x), (b) the directional cosine between the symmetry axis of a fiber and the y-axis (cos y), and (c) the sedimentation velocity (vy) obtained in this study with the results reported by Tian et al. [27]. and Feng et al. [13]. flips from a parallel to a perpendicular orientation, its cross-sectional area decreases, leading to a reduction in the drag force and sudden acceleration. To demonstrate a more robust comparison, the trajectory of the particle in the moving xy plane was plotted and compared (see Figure 9). The particle is driven predominantly in the fluid direction, while gravity governs a slow downward shift. In all cases, the close agreement between our results and those of Tian et al. and Feng et al. demonstrates the correctness of our implementation of the ELER method. The discrepancies may be attributed to differences in numerical setup - the case depends on the viscosity and density of the fluid which are not stated in any of the aforementioned studies. Furthermore, as previously reported by Cui et al. [74], the frequency of flips depends on the Re, which slightly varies between the studies. Different numerical discretization schemes or discrepancies in the lift force formulations can also contribute to variations [13]. 2.4. Grid independence study To assess the validity of the numerical simulations and ensure that the results were not significantly influenced by grid resolution, a grid-independence study was conducted. The grid size was initially chosen based on a previous study [41], which used the same rescaled realistic geometry. To evaluate grid independence, line probes were placed downstream of the first bifurcation near the carina, with two probes in the anterior direction and two in the superior direction (see Figure 11(a)). This area is characterized by turbulent flow conditions, making it suitable for assessing the 17 PREPRINT 0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 4 4.1 4.2 4.3 4.4 4.5·10−4 x(m) y(m) Tian Feng ELER Figure 9: Comparison of the fibers centroid trajectory with the results reported by Tian et al. [27] and Feng et al. [13]. Figure 10: Detail of the computational grid of the right upper lobe (RUL) 18 PREPRINT Table 1: Properties of lattices used for the grid independence study. Number of cells Cell size Time step Relaxation time 40 mil. 0.00018 m 3.5.10−7s 0.5005 60 mil. 0.00015 m 2.6.10−7s 0.5005 82 mil. 0.00013 m 2.0.10−7s 0.5005 impact of grid resolution on flow features. A time period of 0.1 s was simulated with a constant flow rate corresponding to the peak flow rate of the inspiratory cycle (27 l min−1). This flow rate was chosen because the numerical errors tend to be the highest at the maximum fluid velocity. The mean velocity magnitude along the line probes was used for comparison to minimize the influence of local fluctuations. Three different grids with a uniform grid size (for details, see Table 1 and Figure 10) were simulated, and the results are shown in Figure 11(b)-(e). Good agreement in the mean velocity magnitude was observed between the grids with 60 million and 82 million cells. The final simulation was performed using an intermediate grid of 60 million cells. 2.5. Effective diameter simulations The effective diameter simulations are simplified Euler-Lagrange simulations, where the particles are approximated by spheres and the rotation (Equation 11) is neglected. To establish an accurate comparison of ELER with this method, two additional simulations with the most popular empirical models, H-L [22] and T-C [21], on the same geometry under the same breathing conditions were conducted. The hydrodynamic drag force is defined as  FD=1 2CDApρp(u −v)|u −v|(38) where Apis the projected surface area, and CDis the hydrodynamic shape factor. In the H-L correlation, the CDis determined as CD=24 Rep (1 + HaReHb p) + HcRep Hd+ Rep ,(39) where Rep=|u −v|dp ν(40) is the Reynolds number of the particle. νdenotes the fluid kinematic viscosity, and Ha,Hb,Hc and Hdare model-specific parameters that depend on particle sphericity ϕ. T-C defines the shape factor as follows. CD=24 Rep dA deq 1 + 0.15 √cdA deq Rep0.687!+ 0.42 dA deq 2 √c1 + 42500 dA deq Rep−1.16,(41) where dAis the surface equivalent sphere diameter, and cis the surface sphericity. 19 PREPRINT (a) 0 0.2 0.4 0.6 0.8 1 0 1 2 3 4 5 6 Relative length (-) Mean velocity mag. (m/s) 40 mil. 60 mil. 82 mil. (b) 0 0.2 0.4 0.6 0.8 1 0 1 2 3 4 5 6 Relative length (-) Mean velocity mag. (m/s) 40 mil. 60 mil. 82 mil. (c) 0 0.2 0.4 0.6 0.8 1 0 1 2 3 4 5 6 Relative length (-) Mean velocity mag. (m/s) 40 mil. 60 mil. 82 mil. (d) 0 0.2 0.4 0.6 0.8 1 0 1 2 3 4 5 6 Relative length (-) Mean velocity mag. (m/s) 40 mil. 60 mil. 82 mil. (e) Figure 11: Grid independence study: (a) Locations of the line probes used to assess grid independence, (b)-(e) mean velocity magnitude along the line probes for three different grid resolutions; (b) Left branch, horizontal probe (LG1h), (c) Left branch, vertical probe (LG1v), (d) Right branch, horizontal probe (RG1h), (e) Right branch, vertical probe (RG1v). 3. Results and discussion 3.1. Comparison between simulation and experiment Figure 12 compares the numerical and experimental results for the fiber deposition. The simulation data were analyzed at the end of the simulation, and approximately 5% of the particles that had not yet been deposited were excluded from the evaluation. In many segments, good agreement was observed between the simulation and the experiment, particularly in the segments above the bifurcations (segments 1, 2, and 3) and at the funnels at the outlets (segments 13o-22o). Numerical simulation generally predicts a slightly higher deposition fraction in most bifurcations. This discrepancy can be attributed partly to differences in the deposition mechanisms and surface interactions. In the LBM simulations, every contact between a particle and the wall resulted in deposition (the ’perfect sink’ assumption). While the experimental replica surface was coated 20 PREPRINT with silicone oil to promote adhesion and mimic mucus capture [16, 17], the simulation’s assumption of 100% capture upon any contact might still represent an upper bound compared to the in vitro reality, potentially contributing to the observed overestimation. Further investigation of particle-wall interactions under varying surface conditions is necessary to improve the accuracy of the deposition model. Another source of discrepancies may stem from assumptions of linear shear flow in the force calculation. This assumption might result in discrepancies in the areas of the bifurcations ([76]). In addition, the periodic rotational motion of the particles, which is assumed in the ELER method, is not consistently observed in experiments, as reported by Lizal et al. [18]. The largest differences between the simulation and experiment occurred in segments 5, 6, and 7, which corresponded to the second and third generations of branching in the left lung. This excess deposition in the simulation is likely compensated for by a consequent lower deposition fraction in the downstream segments (13, 13o, 16o, 17, 17o), where the experimental results show higher values. The numerical results for the right lung were in good agreement with the experimental data. A further statistics of the experimental data, specifically on counts and dimensions of fibers deposited in each segment can be found in Table A.2 in the Appendix. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 13o14o15o16o17o18o19o20o21o22o 10−3 10−2 10−1 Deposition fraction by count (-) ELER H-L T-C EXP Figure 12: Comparison of deposition fractions obtained from experimental measurements and numerical simulations. To further assess the accuracy of the numerical simulations, a comparison was made between the ELER method and two simplified models based on the effective diameter concept: the T-C model [22] and H-L model [21]. These models approximate the fibers as spheres with equivalent diameters and neglect their rotation. Analysis of the results revealed that the ELER model provided the most accurate predictions of deposition, with the best agreement with the experimental data in 26 segments, compared to five segments for the T-C model and only two segments for the H-L model. The absolute errors for the ELER, T-C, and H-L models were 0.0106, 0.0120, and 0.0146, respectively, confirming the superior accuracy of the ELER approach. This finding is consistent with H-L and T-C comparison studies by Farkas [25]. Based on these results, only the ELER simulation results were considered in the subsequent analysis and discussion. The deposition parameters were analyzed as a function of airway generation to provide a comparison relevant to medical applications. Each segment was assigned its highest generation number, as shown in Figure 2. The deposition fraction and efficiency were then calculated for the following regions: upper airway (UA); trachea (generation 0, G0); generations G1, G2, and G3; and genera21 PREPRINT tions G4-G7. Generations G4-G7 were grouped together because these higher-generation segments are manufactured as single pieces due to their small dimensions. UA G0 G1 G2 G3 G4-G7 10−3 10−2 10−1 100 Deposition efficiency by count (-) ELER EXP Figure 13: Comparison of deposition efficiencies obtained from experimental measurements and numerical simulations, sorted by airway generation. Figure 13 shows the deposition efficiency of each generation. The deposition efficiency generally increased with increasing generation order in both the simulation and experiment, with a steeper increase observed in the LBM results. However, in the experimental results, G2 deviated from this trend, showing a decrease in deposition efficiency compared to G1. Larger discrepancies between the simulation and experiment were observed for the higher generations of branching. The experimental setup used in this study allowed for a unique analysis of the deposition data based on the dimensions of the captured fibers in each segment. This type of analysis is rare in similar experimental studies on fiber deposition. Figure 14 compares the mean values of the volume equivalent diameter (deq) and aspect ratio (β) of the deposited fibers in the upper airways and different generations of branching. In both cases, a reasonable agreement was observed between the simulation and experiment. As expected, particles with higher deq were predominantly deposited in the upper airways and the first few tracheobronchial generations because of their higher inertia. On the other hand, smaller fibers tend to penetrate deeper into the lungs. This trend was evident in both the experimental and numerical results, although the experimental data exhibited slightly more fluctuations. The largest discrepancies between the simulation and experiment occurred in G1, which could be due to the low number of deposited and measured particles in this segment, leading to a higher statistical error. Except for the outlets, the simulations generally predicted a higher mean deq and βthan the experiment. While βincreased with increasing generation downstream of the trachea in both the experimental and numerical results, deq decreased. This observation suggests a complex interplay between the rotational motion of the particles and the tendency of the ELER method to capture particles in bifurcations. However, the flow in the first few bifurcations is turbulent, which makes it challenging to simulate the rotational movements of the fibers accurately. For particles with high β, changes in orientation can significantly alter their trajectory (owing to the dependence of the drag force on the streamwise cross-section, as shown in Equation (13)). This can lead to 22 PREPRINT UA G0 G1 G2 G3 G4-G7 ≧G8 0.25 0.5 0.75 1·10−5 Diameter (m) ELER deq EXP deq 4 8 12 16 Aspect ratio (-) ELER βEXP β Figure 14: Comparison of the mean equivalent diameter (deq ) and mean aspect ratio (β) of deposited fibers for each airway generation. The bar chart (left vertical axis) shows the mean equivalent diameter, while the symbols (right vertical axis) show the mean aspect ratio. differences in the deposition patterns between the simulation and the experiment, particularly for larger particles. This is consistent with the observation that the mean deq of the particles passing through the geometry (i.e., ≥G8) was higher in the experimental setup. A similar behavior was reported by Shachar-Berman et al. [8]. Again, a slight deviation from this trend was observed at G1 in the experiment, likely due to the low number of deposited fibers and the resulting higher statistical error. UA G0 G1 G2 G3 G4-G7 ≧G8 10−3 10−2 10−1 100 Deposition fraction by count (-) ELER 1 µm ELER 2 µm ELER 3 µm EXP 1 µm EXP 2 µm EXP 3 µm Figure 15: Comparison of deposition fractions for fibers sorted into three size groups based on their diameter (dp). A similar trend can be observed in Figure 15, which shows the deposition fraction for the fibers sorted into three classes according to their thickness (dp= 2b= 2c). For small particles (dp= 1 µm and dp= 2 µm), the upper airway and trachea showed good agreement between the 23 PREPRINT simulation and experimental results. However, larger discrepancies were observed in generations G1, G2, and G3, with the differences increasing with the particle thickness. This is consistent with previous observations that orientation changes have a more significant impact on the trajectories of larger particles. In generations G4-G7, the differences between the simulation and experiment were smaller, which could be attributed to the narrower channels in these regions and the transition back to a laminar flow regime. Interestingly, these results differ from those of a previous study using an LBM [41] that simulated spherical particles in a rescaled child airway geometry. In that study, larger discrepancies between the simulation and the experiment were observed for smaller particles. This highlights the impact of the particle shape on the deposition patterns and the importance of considering the fiber orientation in simulations. Figure 16 compares the deposition efficiency obtained in this study with data available in the literature. To facilitate this comparison, the particles were sorted into groups based on their Stokes numbers (Stk) to ensure statistically relevant sample sizes for each evaluation. The trachea and segments up to the fourth generation of branching were evaluated separately. To the best of our knowledge, comparable experimental or simulation data for fiber deposition in a realistic female airway geometry are missing. Therefore, we included data from studies using the male airway geometry by Belka et al. [17] and Farkas et al. [25]. Additionally, we included experimental data from Su et al. [16] and Zhou et al. [15], who used different realistic airway geometries, and Myojo and Takaya [10], who used the idealized Weibel lung model A from the third to fourth generation of branching. Analysis of the data in Figure 16 reveals good agreement between our results and the experimental data for the trachea (Figure 16(a)). Discrepancies compared to the simplified T-C and H-L models are negligible in this region, suggesting less complex flow behavior than in deeper generations. In generations 1-3 (Figures 16(b)–(d)), our simulations predicted a slightly higher deposition efficiency than the experimental data. This difference may be due to the use of a realistic inspiration profile in our simulations, which includes a higher peak velocity that could increase deposition in comparison with studies that used constant flowrates. The differences between the ELER model and simplified models become more pronounced in higher generations, highlighting the importance of considering fiber rotation in the branching regions. It is worth noting that the experimental data from Su et al. [16] have a higher standard deviation, which should be considered when comparing our results. Overall, the ELER simulations tend to overestimate the fiber deposition compared with our experiments, which is consistent with previous findings. In Generation 4 (Figure 16(e)), both our experimental and numerical results show higher deposition across the entire Stk range compared to Myojo and Takaya [10]. This discrepancy could be attributed to the use of an idealized geometry and the absence of upstream bifurcations in their model, which does not account for the full flow history. The experiments of Zhou et al. [15] were conducted for higher Stk values; however, extrapolating their data to the lower Stk range would still indicate lower deposition than that observed in our study. Interestingly, the experimental data from Belka et al. [17] show a lower deposition efficiency for almost all segments and Stk numbers compared to other studies, including the corresponding H-L and T-C simulations performed by Farkas et al. [25]. This discrepancy could be related to differences in the experimental setup or specific characteristics of the airway models used. In summary, the female airway geometry used in this study appears to result in a slightly higher deposition rate compared to the male airway geometry, with a shift of the capture-efficiency curve 24 PREPRINT •The initial position and orientation of the particles released in the simulations were randomly generated, as it is challenging to measure these parameters experimentally. This uncertainty in the initial conditions could have affected the deposition patterns and introduced inaccuracies in the results. •The subgrid-scale turbulence fluctuations smaller than the grid size were not considered in this study and can lead to small discrepancies in trajectories and orientations, and can influence exact deposition positions for particles near the walls, predominantly when reaching the peak flow rate. Modeling of the turbulence fluctuations, using, e.g., a Gaussian random filter ([56]), will be involved in future work investigating higher flow rates. •A dilute air and one-way coupling model was used, assuming that the particles do not influence the airflow. This assumption may not be valid in all regions of the airway, particularly in the turbulent flow regime near bifurcations, where two-way coupling might be necessary for higher accuracy. •Particle concentration was assumed to be low, and inter-particle collisions were neglected. Additionally, potential changes in the particle shape owing to collisions were not considered. •Due to the time-consuming nature of the experiments, only one repetition was performed. Additional experimental data would be useful for a better estimation of the uncertainties. •The effect of expiration on particle deposition was not considered in this study. This could be particularly relevant for particles deposited towards the end of the inspiration cycle, as they might be exhaled before having a chance to deposit permanently. •Only one flow rate was simulated in this study. Systematic exploration of the influence of various breathing patterns was performed by e.g. [80, 81], and hence similar trends in the deposition characteristics can be expected. 5. Conclusion This study successfully applied the LBM in conjunction with the ELER method to simulate fiber transport and deposition in a realistic female airway model. A unique comparison of the deposition characteristics between numerical simulations and experiments conducted on the same airway geometry, extending to the 7th generation of branching, is presented. The results showed good agreement in the upper airways and trachea, but some discrepancies were observed in the bronchial bifurcations, likely owing to the challenges in capturing complex flow phenomena (such as Dean vortices and turbulence) using the ELER method. In addition to comparing the ELER method with experimental data, its accuracy was also assessed by comparing it with simplified fiber deposition models that neglect fiber rotation (T-C and H-L methods). The ELER outperformed these simplified models, demonstrating its superior ability to predict fiber deposition in a realistic airway geometry and emphasizing the importance of orientation-dependent calculation of fiber transport and deposition. This study also enabled a detailed analysis of the deposition patterns for different fiber size groups. Notably, the ELER method exhibited better agreement with the experimental data for smaller particles, indicating a higher accuracy in predicting the deposition fractions for these size ranges. Furthermore, the results confirmed that fibers with higher aspect ratios had a greater tendency to penetrate deeper into the lungs, which is consistent with previous findings. A comparison of the deposition efficiency with literature data, using Stk as a basis for comparison, showed good agreement in the trachea. However, discrepancies were observed in the bifurcations, potentially owing to the differences in airway geometries and flow regimes between 31 PREPRINT the studies. The use of a realistic, transient inspiration profile in our simulations led to a higher deposition efficiency compared to studies that employed steady flow conditions, likely owing to the higher peak velocities and resulting higher Stk numbers in our simulations. The numerical analysis revealed that particles released at the beginning of the inspiration cycle, when the flow rates were the highest, were more likely to deposit in the upper airways. This finding highlights the importance of considering initial high flow rates when designing drug delivery strategies. For instance, adjusting the release timing to avoid this initial phase could help prevent unwanted deposition in the upper airways. No significant differences in deposition were observed for particles released after the first third of the inspiration cycle, suggesting that the initial high flow rates had the greatest impact on the deposition location. Furthermore, only negligible differences in deposition were found among the various branching segments, indicating a relatively uniform distribution of deposition throughout the tracheobronchial tree. Approximately 70% of the particles were deposited during the first inspiration cycle. Notably, the orientation-dependent deposition mechanism, which accounts for the angle of the fiber relative to the airway wall, was activated in 93% of the particles. This emphasizes its crucial role in accurately predicting deposition, particularly for applications such as locally targeted drug delivery, in which precise control over the deposition location is essential. Simulating fibrous particles in realistic airway geometries using a realistic breathing profile provides valuable insights into optimizing the inhalation process and predicting fiber deposition patterns. However, the results also demonstrated that current numerical techniques for fiber modeling, specifically ELER, and the drag models H-L and T-C, still exhibit some discrepancies compared to the experimental data, especially in the branching regions. Further research is needed to improve these techniques and expand their applicability, ultimately enabling them to replace or complement experimental studies in a wider range of scenarios and contributing to the ongoing development of effective medical treatments. This study also highlights a key difference between spherical and fibrous particles: fibers exhibit a greater ability to bypass the upper airways and penetrate deeper into the lungs. This finding has important implications for understanding the potential health risks associated with inhalation of different types of particles and for developing targeted drug delivery strategies that can effectively reach specific regions of the respiratory system. Finally, it is crucial to consider the influence of individual airway geometries and breathing patterns on particle deposition. This study focused on a female airway model, but future research should investigate the differences in deposition between male, female, and child airways and evaluate the impact of various inhalation regimes. This knowledge is essential for developing personalized inhalation therapies tailored to the specific characteristics of individual patients. Acknowledgment This work was supported by the Ministry of Education, Youth and Sports of the Czech Republic through the e-INFRA CZ (ID: 90254), by the Czech Science Foundation grant 22-20357S, and by the internal research project of Brno University of Technology, reg. no. FSI-S-23-8192. The authors express gratitude to Pavel Fr¨oml for conducting a series of experimental measurements. 32 PREPRINT Declaration of generative AI and AI-assisted technologies in the writing process During the preparation of this study, the authors used ChatGPT, Grammarly, and PaperPal to improve the language and readability. After using these tools, the authors reviewed and edited the content as needed and took full responsibility for the contents of the published article. Data availability The data used in this article are accessible in the Zenodo repository under the following DOI: 10.5281/zenodo.15166271. 33 PREPRINT Appendix A. Detail to the experimental setup Figure A.22: A representative micrograph showing collected glass fibers on a filter membrane 34 PREPRINT Table A.2: Summary statistics of the total fiber counts and their length 2aand thickness dpobtained by the experimental measurement. Segment no. No. of particles Mean of dp(µm) Std. deviation of dp(µm) Mean of 2a(µm) Std. deviation of 2a(µm) 1 325264 2.0 1.3 15.8 14.5 2 96341 2.4 1.2 18.2 15.3 3 37160 1.9 1.0 20.1 15.9 4 42206 1.7 0.9 14.7 8.3 5 43583 1.9 1.1 14.8 11.5 6 42106 2.2 1.8 17.0 14.0 7 36242 1.6 0.6 13.7 9.4 8 29820 2.4 1.8 13.9 9.6 9 26608 1.4 0.6 14.3 8.3 10 20644 1.7 0.9 13.7 11.6 11 18351 1.6 0.7 21.3 18.5 12 22021 2.2 1.9 21.7 15.9 13 83954 1.9 0.7 14.3 10.3 14 30278 1.4 0.6 21.8 28.9 15 22021 1.9 1.1 15.7 12.8 16 112397 1.8 0.6 15.9 11.7 17 26608 1.9 0.8 21.2 13.7 18 53676 1.8 0.9 19.2 24.6 19 61016 1.9 0.7 21.4 24.4 20 23397 1.3 0.5 13.6 7.6 21 37160 1.5 0.8 15.1 10.6 22 38995 1.7 0.5 13.8 9.6 13o 856273 1.5 0.6 17.4 15.9 14o 585743 1.5 0.7 16.2 13.9 15o 757403 1.6 0.7 16.7 14.8 16o 1031666 1.6 0.7 17.3 15.2 17o 340679 1.5 0.6 16.5 16.9 18o 1087424 1.6 0.8 17.4 16.3 19o 837929 1.6 0.7 16.3 14.3 20o 347844 1.5 0.7 16.8 15.0 21o 847405 1.7 0.7 16.8 15.6 22o 1047117 1.6 0.7 16.1 14.9 References [1] O. P. Soldin, D. R. Mattison, Sex differences in pharmacokinetics and pharmacodynamics, Clinical Pharmacokinetics 48 (3) (2009) 143 – 157. doi:10.2165/00003088-200948030-00001. [2] C. M. Madla, F. K. Gavins, H. A. Merchant, M. Orlu, S. Murdan, A. W. Basit, Let’s talk about sex: Differences in drug therapy in males and females, Advanced Drug Delivery Reviews 175 (2021) 113804. doi:10.1016/j. addr.2021.05.014. [3] S. Sharifi, G. Caracciolo, D. Pozzi, L. Digiacomo, J. Swann, H. E. Daldrup-Link, M. Mahmoudi, The role of sex as a biological variable in the efficacy and toxicity of therapeutic nanomedicine, Advanced Drug Delivery Reviews 174 (2021) 337–347. doi:oi.org/10.1016/j.addr.2021.04.028. [4] S. Cheng, J. Butler, S. Gandevia, L. Bilston, Movement of the tongue during normal breathing in awake healthy humans, Journal of Physiology 586 (17) (2008) 4283 – 4294. doi:10.1113/jphysiol.2008.156430. [5] A. LoMauro, A. Aliverti, Sex differences in respiratory function, Breathe 14 (2) (2018) 131–140. doi:10.1183/ 20734735.000318. [6] M. A. Carey, J. W. Card, J. W. Voltz, D. R. Germolec, K. S. Korach, D. C. Zeldin, The impact of sex and sex hormones on lung physiology and disease: lessons from animal studies, American Journal of Physiology-Lung Cellular and Molecular Physiology 293 (2) (2007) L272–L278. doi:10.1152/ajplung.00174.2007. [7] T. Mekonnen, X. Cai, C. Burchell, H. Gholizadeh, S. Cheng, A review of upper airway physiology relevant to the delivery and deposition of inhalation aerosols, Advanced Drug Delivery Reviews 191 (2022) 114530. doi:10.1016/j.addr.2022.114530. 35 PREPRINT [8] L. Shachar-Berman, S. Bhardwaj, Y. Ostrovski, P. Das, P. Koullapis, S. Kassinos, J. Sznitman, In Silico Optimization of Fiber-Shaped Aerosols in Inhalation Therapy for Augmented Targeting and Deposition across the Respiratory Tract, Pharmaceutics 12 (3) (2020) 230. doi:10.3390/pharmaceutics12030230. [9] J. Marijnlssen, A. Zeckendorf, S. Lemkowltz, H. Bibo, Transport and deposition of uniform respirable fibres in a physical lung model, Journal of Aerosol Science 22 (1991) S859–S862. doi:10.1016/S0021-8502(05)80234-4. [10] T. Myojo, M. Takaya, Estimation of fibrous aerosol deposition in upper bronchi based on experimental data with model bifurcation, Industrial health 39 (2) (2001) 141–149. doi:10.2486/indhealth.39.141. [11] E. R. Weibel, A. F. Cournand, D. W. Richards, Morphometry of the human lung, Vol. 1, Springer, 1963. doi:10.1007/978-3-642-87553-3. [12] W.-C. Su, Y. S. Cheng, Deposition of fiber in a human airway replica, Journal of aerosol science 37 (11) (2006) 1429–1441. doi:j.jaerosci.2006.01.015. [13] Y. Feng, C. Kleinstreuer, Analysis of Non-Spherical Particle Transport in Complex Internal Shear Flows, Physics of Fluids 25 (2013) 1904. doi:10.1063/1.4821812. [14] W.-C. Su, Y. S. Cheng, Fiber deposition pattern in two human respiratory tract replicas, Inhalation toxicology 18 (10) (2006) 749–760. doi:10.1080/08958370600748513. [15] Y. Zhou, W.-C. Su, Y. S. Cheng, Fiber deposition in the tracheobronchial region: Experimental measurements, Inhalation toxicology 19 (13) (2007) 1071–1078. doi:10.1080/08958370701626634. [16] W.-C. Su, Y. S. Cheng, Deposition of man-made fibers in human respiratory airway casts, Journal of aerosol science 40 (3) (2009) 270–284. doi:10.1016/j.jaerosci.2008.11.003. [17] M. Belka, F. Lizal, J. Jedelsky, J. Elcner, P. K. Hopke, M. Jicha, Deposition of glass fibers in a physically realistic replica of the human respiratory tract, Journal of Aerosol Science 117 (2018) 149–163. doi:10.1016/ j.jaerosci.2017.11.006. [18] F. Lizal, M. Cabalka, M. Maly, J. Elcner, M. Belka, E. Lizalova Sujanska, A. Farkas, P. Starha, O. Pech, O. Misik, J. Jedelsky, M. Jicha, On the behavior of inhaled fibers in a replica of the first airway bifurcation under steady flow conditions, Aerosol Science and Technology 56 (4) (2022) 367–381. doi:10.1080/02786826.2022.2027334. [19] C. Kleinstreuer, Y. Feng, Computational Analysis of Non-Spherical Particle Transport and Deposition in Shear Flow With Application to Lung Aerosol Dynamics—A Review, Journal of Biomechanical Engineering 135 (021008) (Feb. 2013). doi:10.1115/1.4023236. [20] I. A. Lasso, P. Weidman, Stokes drag on hollow cylinders and conglomerates, The Physics of fluids 29 (12) (1986) 3921–3934. doi:10.1063/1.865732. [21] A. Haider, O. Levenspiel, Drag coefficient and terminal velocity of spherical and nonspherical particles, Powder technology 58 (1) (1989) 63–70. doi:10.1016/0032-5910(89)80008-7. [22] S. Tran-Cong, M. Gay, E. E. Michaelides, Drag coefficients of irregularly shaped particles, Powder Technology 139 (1) (2004) 21–32. doi:10.1016/j.powtec.2003.10.002. [23] W. St¨ober, Dynamic shape factors of nonspherical aerosol particles, Assessment of airborne particles (1972) 249–289. [24] K. Inthavong, J. Wen, Z. Tian, J. Tu, Numerical study of fibre deposition in a human nasal cavity, Journal of Aerosol Science 39 (2008) 253–265. doi:10.1016/j.jaerosci.2007.11.007. [25] A. Farkas, F. Lizal, J. Elcner, J. Jedelsky, M. Jicha, Numerical simulation of fibre deposition in oral and large bronchial airways in comparison with experiments, Journal of Aerosol Science 136 (2019) 1–14. doi: 10.1016/j.jaerosci.2019.06.003. [26] X. Chen, W. Zhong, J. Tom, C. Kleinstreuer, Y. Feng, X. He, Experimental-computational study of fibrous particle transport and deposition in a bifurcating lung model, Particuology 28 (2016) 102–113. doi:10.1016/ j.partic.2016.02.002. [27] L. Tian, G. Ahmadi, Z. Wang, P. K. Hopke, Transport and deposition of ellipsoidal fibers in low Reynolds number flows, Journal of Aerosol Science 45 (2012) 1–18. doi:10.1016/j.jaerosci.2011.09.001. [28] L. Tian, G. Ahmadi, Fiber transport and deposition in human upper tracheobronchial airways, Journal of Aerosol Science (60) (2013) 1–20. doi:10.1016/j.jaerosci.2013.02.001. [29] K. T. Shanley, G. Ahmadi, P. K. Hopke, Y.-S. Cheng, Simulated airflow and rigid fiber behavior in a realistic nasal airway model, Particulate Science and Technology 36 (2) (2018) 131–140, 10.1080/02726351.2016.1208694. doi:10.1080/02726351.2016.1208694. [30] J. Li, J. Ma, J. Dong, W. Yang, G. Ahmadi, J. Tu, L. Tian, Microfiber transport characterization in human nasal cavity – Effect of fiber length, Journal of Aerosol Science 160 (2022) 105908. doi:10.1016/j.jaerosci. 2021.105908. [31] J. Li, J. Ma, G. Ahmadi, J. Dong, W. Yang, J. Tu, L. Tian, Shear induced lift and rotation on MicroFiber deposition in low Reynolds number flows, Journal of Aerosol Science 167 (2023) 106094. doi:10.1016/j. 36 PREPRINT jaerosci.2022.106094. [32] J. Li, J. Ma, J. Dong, W. Yang, J. Tu, L. Tian, Total and regional microfiber transport characterization in a 15th - Generation human respiratory airway, Computers in Biology and Medicine 163 (2023) 107180. doi:10.1016/j.compbiomed.2023.107180. [33] M. Kiasadegh, H. Emdad, G. Ahmadi, O. Abouali, Transient numerical simulation of airflow and fibrous particles in a human upper airway model, Journal of Aerosol Science 140 (2020) 105480. doi:10.1016/j.jaerosci.2019. 105480. [34] M. M. Tavakol, E. Ghahramani, O. Abouali, M. Yaghoubi, G. Ahmadi, Deposition fraction of ellipsoidal fibers in a model of human nasal cavity for laminar and turbulent flows, Journal of Aerosol Science 113 (2017) 52–70. doi:10.1016/j.jaerosci.2017.07.008. [35] M. Abolhassantash, M. M. Tavakol, O. Abouali, M. Yaghoubi, G. Ahmadi, Deposition fraction of ellipsoidal fibers in the human nasal cavityInfluence of non-creeping formulation of hydrodynamic forces and torques, International Journal of Multiphase Flow 126 (2020) 103238. doi:10.1016/j.ijmultiphaseflow.2020.103238. [36] M. Zastawny, G. Mallouppas, F. Zhao, B. van Wachem, Derivation of drag and lift force and torque coefficients for non-spherical particles in flows, International Journal of Multiphase Flow 39 (2012) 227–239. doi:10.1016/ j.ijmultiphaseflow.2011.09.004. [37] R. Ouchene, M. Khalij, B. Arcen, A. Tani`ere, A new set of correlations of drag, lift and torque coefficients for non-spherical particles and large Reynolds numbers, Powder Technology 303 (2016) 33–43. doi:10.1016/j. powtec.2016.07.067. [38] A. A. Mofakham, G. Ahmadi, On random walk models for simulation of particle-laden turbulent flows, International Journal of Multiphase Flow 122 (2020) 103157. doi:10.1016/j.ijmultiphaseflow.2019.103157. [39] L. Tian, G. Ahmadi, Computational modeling of fiber transport in human respiratory airways—A review, Experimental and Computational Multiphase Flow 3 (1) (2021) 1–20. doi:10.1007/s42757-020-0061-7. [40] T. Henn, G. Th¨ater, W. D¨orfler, H. Nirschl, M. Krause, Parallel dilute particulate flow simulations in the human nasal cavity, Computers & Fluids 124 (2016) 197–207. doi:10.1016/j.compfluid.2015.08.002. [41] F. Prinz, J. Pokorn´y, J. Elcner, F. L´ızal, O. Miˇs´ık, M. Mal´y, M. Bˇelka, N. Hafen, A. Kummerl¨ander, M. J. Krause, J. Jedelsk´y, M. J´ıcha, Comprehensive experimental and numerical validation of Lattice Boltzmann fluid flow and particle simulations in a child respiratory tract, Computers in Biology and Medicine 170 (2024) 107994. doi:10.1016/j.compbiomed.2024.107994. [42] F. Lizal, J. Elcner, P. K. Hopke, J. Jedelsky, M. Jicha, Development of a realistic human airway model, Proceedings of the Institution of Mechanical Engineers, Part H: Journal of Engineering in Medicine 226 (3) (2012) 197–207. doi:10.1177/0954411911430188. [43] N. Griscom, M. Wohl, Dimensions of the growing trachea related to age and gender, American Journal of Roentgenology 146 (2) (1986) 233–237, pMID: 3484568. doi:10.2214/ajr.146.2.233. [44] T. R. Martin, R. G. Castile, J. J. Fredberg, M. E. Wohl, J. Mead, Airway size is related to sex but not lung size in normal adults, Journal of Applied Physiology 63 (5) (1987) 2042–2047, publisher: American Physiological Society. doi:10.1152/jappl.1987.63.5.2042. [45] A. W. Sheel, J. A. Guenette, R. Yuan, L. Holy, J. R. Mayo, A. M. McWilliams, S. Lam, H. O. Coxson, Evidence for dysanapsis using computed tomographic imaging of the airways in older ex-smokers, Journal of Applied Physiology 107 (5) (2009) 1622–1628, publisher: American Physiological Society. doi:10.1152/japplphysiol. 00562.2009. [46] ICRP, Human respiratory tract model for radiological protection, Annals of the ICRP 66 (1994) 1–3. [47] N. Jahani, S. Choi, J. Choi, K. Iyer, E. A. Hoffman, C.-L. Lin, Assessment of regional ventilation and deformation using 4d-ct imaging for healthy human lungs during tidal breathing, Journal of Applied Physiology 119 (10) (2015) 1064–1074. doi:10.1152/japplphysiol.00339.2015. [48] W. H. Organization, et al., Determination of airborne fibre number concentrations: a recommended method, by phase-contrast optical microscopy (membrane filter method), World Health Organization, 1997. [49] T. Kr¨uger, H. Kusumaatmaja, A. Kuzmin, O. Shardt, G. Silva, E. M. Viggen, The Lattice Boltzmann Method: Principles and Practice, Graduate Texts in Physics, Springer International Publishing. doi: 10.1007/978-3-319-44649-3. [50] A. Kummerl¨ander, F. Bukreev, S. F. R. Berg, M. Dorn, M. J. Krause, Advances in computational process engineering using lattice boltzmann methods on high performance computers, in: W. E. Nagel, D. H. Kr¨oner, M. M. Resch (Eds.), High Performance Computing in Science and Engineering ’22, Springer Nature Switzerland, pp. 233–247. doi:10.1007/978-3-031-46870-4\_16. [51] M. Haussmann, F. Ries, J. B. Jeppener-Haltenhoff, Y. Li, M. Schmidt, C. Welch, L. Illmann, B. B¨ohm, H. Nirschl, M. J. Krause, A. Sadiki, Evaluation of a near-wall-modeled large eddy lattice boltzmann method for 37 PREPRINT the analysis of complex flows relevant to IC engines 8 (2) 43. doi:10.3390/computation8020043. [52] S. Chapman, T. G. Cowling, The mathematical theory of non-uniform gases: an account of the kinetic theory of viscosity, thermal conduction and diffusion in gases, Cambridge university press, 1990. [53] S. Hou, J. D. Sterling, S. Chen, G. D. Doolen, A lattice boltzmann subgrid model for high reynolds number flows, arXiv: Cellular Automata and Lattice Gases (1994). [54] I. B. Celik, Z. N. Cehreli, I. Yavuz, Index of Resolution Quality for Large Eddy Simulations, Journal of Fluids Engineering 127 (5) (2005) 949–958. doi:10.1115/1.1990201. URL https://doi.org/10.1115/1.1990201 [55] M. Sommerfeld, O. L. Sgrott, M. A. Taborda, P. Koullapis, K. Bauer, S. Kassinos, Analysis of flow field and turbulence predictions in a lung model applying RANS and implications for particle deposition, European Journal of Pharmaceutical Sciences 166 (2021) 105959. doi:10.1016/j.ejps.2021.105959. URL https://www.sciencedirect.com/science/article/pii/S0928098721002621 [56] M. Salmanzadeh, M. Rahnama, G. Ahmadi, Effect of Sub-Grid Scales on Large Eddy Simulation of Particle Deposition in a Turbulent Channel Flow, Aerosol Science and Technology 44 (9) (2010) 796–806. doi:10.1080/ 02786826.2010.492052. [57] H. Sajjadi, M. Salmanzadeh, G. Ahmadi, S. Jafari, Simulations of indoor airflow and particle dispersion and deposition by the lattice Boltzmann method using LES and RANS approaches, Building and Environment 102 (2016) 1–12. doi:https://doi.org/10.1016/j.buildenv.2016.03.006. URL https://www.sciencedirect.com/science/article/pii/S0360132316300762 [58] E. Ghahramani, O. Abouali, H. Emdad, G. Ahmadi, Numerical investigation of turbulent airflow and microparticle deposition in a realistic model of human upper airway using LES, Computers & Fluids 157 (2017) 43–54. doi:10.1016/j.compfluid.2017.08.003. URL https://www.sciencedirect.com/science/article/pii/S0045793017302748 [59] P. Koullapis, S. C. Kassinos, J. Muela, C. Perez-Segarra, J. Rigola, O. Lehmkuhl, Y. Cui, M. Sommerfeld, J. Elcner, M. Jicha, I. Saveljic, N. Filipovic, F. Lizal, L. Nicolaou, Regional aerosol deposition in the human airways: The SimInhale benchmark case and a critical assessment of in silico methods, European Journal of Pharmaceutical Sciences 113 (2018) 77–94. doi:10.1016/j.ejps.2017.09.003. URL https://www.sciencedirect.com/science/article/pii/S0928098717304992 [60] P. Skordos, Initial and boundary conditions for the lattice Boltzmann method, Phys. Rev. E 48 (1993) 4823 – 4842. doi:10.1103/PhysRevE.48.4823. [61] M. Bouzidi, M. Firdaouss, P. Lallemand, Momentum transfer of a Boltzmann-lattice fluid with boundaries, Phys. Fluids 13 (2001) 3452–3459. doi:10.1063/1.1399290. [62] A. Lintermann, W. Schr¨oder, Simulation of aerosol particle deposition in the upper human tracheobronchial tract, European Journal of Mechanics - B/Fluids 63 (2017) 73–89. doi:10.1016/j.euromechflu.2017.01.008. [63] F.-G. Fan, G. Ahmadi, A sublayer model for wall deposition of ellipsoidal particles in turbulent streams, Journal of Aerosol Science 26 (5) (1995) 813–840. doi:10.1016/0021-8502(95)00021-4. [64] H. Brenner, The Stokes resistance of an arbitrary particle—IV Arbitrary fields of flow, Chemical Engineering Science 19 (10) (1964) 703–727. doi:10.1016/0009-2509(64)85084-3. [65] E. Y. Harper, I.-D. Chang, Maximum dissipation resulting from lift in a slow viscous shear flow, Journal of Fluid Mechanics 33 (2) (1968) 209–225. doi:10.1017/S0022112068001254. [66] Y. Cui, J. Ravnik, M. Hriberˇsek, P. Steinmann, On Constitutive Models for the Momentum Transfer to Particles in Fluid-Dominated Two-Phase Flows, in: H. Altenbach, F. Jablonski, W. H. M¨uller, K. Naumenko, P. Schneider (Eds.), Advances in Mechanics of Materials and Structural Analysis: In Honor of Reinhold Kienzler, Advanced Structured Materials, Springer International Publishing, Cham, 2018, pp. 1–25. doi:10.1007/978-3-319-70563-7\_1. [67] Y. Cui, J. Ravnik, M. Hriberˇsek, P. Steinmann, Towards a unified shear-induced lift model for prolate spheroidal particles moving in arbitrary non-uniform flow, Computers & Fluids 196 (2020) 104323. doi:10.1016/j. compfluid.2019.104323. [68] G. B. Jeffery, The motion of ellipsoidal particles immersed in a viscous fluid, Proceedings of the Royal Society of London. Series A, Containing papers of a mathematical and physical character 102 (715) (1922) 161–179. [69] F.-G. Fan, G. Ahmadi, A sublayer model for wall deposition of ellipsoidal particles in turbulent streams, Journal of Aerosol Science 26 (5) (1995) 813–840. doi:10.1016/0021-8502(95)00021-4. [70] A. Kummerlander, T. Bingert, F. Bukreev, L. E. Czelusniak, D. Dapelo, N. Hafen, M. Heinzelmann, S. Ito, J. Jeßberger, H. Kusumaatmaja, J. E. Marquardt, M. Rennick, T. Pertzel, F. Prinz, M. Sadric, M. Schecher, S. Simonis, P. Sitter, D. Teutscher, M. Zhong, M. J. Krause, OpenLB Release 1.7: Open Source Lattice Boltzmann Code (Feb. 2024). doi:10.5281/zenodo.10684609. 38 PREPRINT [71] A. Kummerlander, T. Bingert, F. Bukreev, L. E. Czelusniak, D. Dapelo, S. Englert, N. Hafen, M. Heinzelmann, S. Ito, J. Je˜ A¨ Yberger, F. Kaiser, E. Kummer, H. Kusumaatmaja, J. E. Marquardt, M. Rennick, T. Pertzel, F. Prinz, M. Sadric, M. Schecher, S. Simonis, P. Sitter, D. Teutscher, M. Zhong, M. J. Krause, Openlb user guide 1.7 (Aug. 2024). doi:10.5281/zenodo.13293033. [72] M. Krause, A. Kummerl¨ander, S. Avis, H. Kusumaatmaja, D. Dapelo, F. Klemens, M. Gaedtke, N. Hafen, A. Mink, R. Trunk, J. Marquardt, M. Maier, M. Haussmann, S. Simonis, OpenLB–Open source lattice Boltzmann code, Computers & Mathematics with Applications 81 (2021) 258–288. doi:10.1016/j.camwa.2020.04. 033. [73] K. Shanley, G. Ahmadi, A Numerical Model for Simulating the Motions of Ellipsoidal Fibers Suspended in Low Reynolds Number Shear Flows, Aerosol Science and Technology 45 (2011) 838–848. doi:10.1080/02786826. 2011.566293. [74] Y. Cui, J. Ravnik, M. Hriberˇsek, P. Steinmann, A novel model for the lift force acting on a prolate spheroidal particle in an arbitrary non-uniform flow. Part I. Lift force due to the streamwise flow shear, International Journal of Multiphase Flow 104 (2018) 103–112. doi:10.1016/j.ijmultiphaseflow.2018.03.007. [75] J. Wedel, P. Steinmann, M. ˇ Strakl, M. Hriberˇsek, J. Ravnik, Shape matters: Lagrangian tracking of complex nonspherical microparticles in superellipsoidal approximation, International Journal of Multiphase Flow 158 (2023) 104283. doi:10.1016/j.ijmultiphaseflow.2022.104283. [76] J. Elcner, F. Lizal, J. Jedelsky, M. Jicha, M. Chovancova, Numerical investigation of inspiratory airflow in a realistic model of the human tracheobronchial airways and a comparison with experimental results, Biomechanics and Modeling in Mechanobiology 15 (2) (2016) 447–469. doi:10.1007/s10237-015-0701-1. [77] A. Naseri, S. Shaghaghian, O. Abouali, G. Ahmadi, Numerical investigation of transient transport and deposition of microparticles under unsteady inspiratory flow in human upper airways, Respiratory Physiology & Neurobiology 244 (2017) 56–72. doi:10.1016/j.resp.2017.06.005. [78] H. Bahmanzadeh, O. Abouali, G. Ahmadi, Unsteady particle tracking of micro-particle deposition in the human nasal cavity under cyclic inspiratory flow, Journal of Aerosol Science 101 (2016) 86–103. doi:10.1016/j. jaerosci.2016.07.010. [79] J. Wedel, P. Steinmann, F. Prinz, F. L´ızal, M. Hriberˇsek, J. Ravnik, Mass distribution impacts on particle translation and orientation dynamics in dilute flows, Powder Technology 452 (2025) 120424. doi:10.1016/j. powtec.2024.120424. [80] K. Kuga, R. Kizuka, N. D. Khoa, K. Ito, Effect of transient breathing cycle on the deposition of micro and nanoparticles on respiratory walls, Computer Methods and Programs in Biomedicine 236 (2023) 107501. doi: 10.1016/j.cmpb.2023.107501. [81] H. Liu, S. Ma, T. Hu, D. Ma, Computational investigation of flow characteristics and particle deposition patterns in a realistic human airway model under different breathing conditions, Respiratory Physiology & Neurobiology 314 (2023) 104085. doi:10.1016/j.resp.2023.104085. 39