key: cord-009766-mdmqcvww
authors: Comerford, Andrew; Förster, Christiane; Wall, Wolfgang A.
title: Structured Tree Impedance Outflow Boundary Conditions for 3D Lung Simulations
date: 2010-08-01
journal: J Biomech Eng
DOI: 10.1115/1.4001679
sha: 
doc_id: 9766
cord_uid: mdmqcvww

In this paper, we develop structured tree outflow boundary conditions for modeling the airflow in patient specific human lungs. The utilized structured tree is used to represent the nonimageable vessels beyond the 3D domain. The coupling of the two different scales (1D and 3D) employs a Dirichlet–Neumann approach. The simulations are performed under a variety of conditions such as light breathing and constant flow ventilation (which is characterized by very rapid acceleration and deceleration). All results show that the peripheral vessels significantly impact the pressure, however, the flow is relatively unaffected, reinforcing the fact that the majority of the lung impedance is due to the lower generations rather than the peripheral vessels. Furthermore, simulations of a hypothetical diseased lung (restricted flow in the superior left lobe) under mechanical ventilation show that the mean pressure at the outlets of the 3D domain is about 28% higher. This hypothetical model illustrates potential causes of volutrauma in the human lung and furthermore demonstrates how different clinical scenarios can be studied without the need to assume the unknown flow distribution into the downstream region.

Mechanical ventilation is an indispensable tool for the survival of critical care patients suffering from acute respiratory distress syndrome ͑ARDS͒ or acute lung injury ͑ALI͒. However, the mortality rates of these diseases remains high, approximately 40% ͓1͔. In addition, it is widely known that the use of mechanical ventilation is itself the cause of a number of further associated complications ͓2͔, which are collectively termed ventilator induced lung injury ͑VILI͒. This is despite the more recent adoption of modern protective ventilation strategies, where lower tidal volumes are administered while retaining a positive pressure at the end of expiration, which increases lung functional residual capacity and alveolar recruitment. These protocols improve the situation considerably over conventional ventilation ͓3-5͔. However, there are competing ventilator strategies each with their own merits and drawbacks, for example, the level of positive end-expiratory pressure ͑PEEP͒ required by the patient ͓6,7͔. Biologically, these aforementioned ventilator related diseases manifest themselves at the alveolar level and are characterized by inflammation of the lung parenchyma, albeit still much controversy surrounds the precise mechanisms ͓8͔. VILI can be characterized into separate but not mutually exclusive protocol dependent mechanisms: these are volutrauma ͑alveolar overdistension͒, atelectotrauma ͑recruitment-derecruitment͒, and biotrauma ͑inflammation͒ ͓9͔. Predisposition to these conditions is brought about by heterogeneous lung mechanics when the patient is suffering from ARDS or ALI ͓9͔, for example, the functional capacity of a diseased lung can be reduced resulting in alveolar overdistension ͓10͔. Mortality is usually caused via indirect mechanisms resulting from VILI ͓11͔ in which a cascade of events leads to sepsis or multiple organ failure. Understanding the reason why the lungs still become damaged or inflamed during mechanical ventilation is a key question sought by the medical community.

Computational modeling of the tracheobronchial region repre-sents an approach that can provide more in depth knowledge into how the human lung functions mechanically, both during normal breathing and with the assistance of mechanical ventilation. Of particular interest are what conditions within the lung mechanically differ under both diseased and healthy conditions. Currently, a number of limitations exist for modeling the lower airways ͑in this paper lower airways and higher generations are used synonymously to mean peripheral vessels͒. First the resolution of CT imaging is restricted ͑maximum 0.5 mm resolution͒, hence, smaller vessels are not visible for segmentation. Second, even if the entire lung was fully segmentable, currently, it is not computationally feasible to simulate the full lung tree. Due to these points, alternative methods for modeling the lung must be sought. Previously, there have been extensive models of nonpatient specific human lungs ͓12-14͔ based around literature data such as Weibel ͓15͔ and Horsfield et al. ͓16͔ . Furthermore, more recent studies have utilized patient specific CT data ͓17-19͔ and CT ovine lungs ͓20͔ under a variety of boundary conditions, including zero pressure conditions, prescribed flowrate, and impedance. In addition, the simulations of Refs. ͓17,21͔ included for the first time the effects of fluid-structure interaction ͑FSI͒.

There has been considerable work in the area of lung impedance, specifically related to 1D transmission line models ͓22-25͔. Impedance of the lung is an extremely important phenomenon as it plays a central role in the development and distribution of lung disease. More specifically, the heterogeneity of lung disease is brought about by changes in impedance conditions within the lung. These changes further affect how flow and pressure distribute through the lung. The aforementioned impedance models assume 1D fluid mechanics throughout the entire lung. However, from previous 3D simulations it is well known that the flow in the upper portion of the lower airways is not accurately modeled with 1D equations since typical bifurcating flow phenomena, such as flow recirculation, exist. For this reason, a fully coupled 1D-3D approach must be used and is hence developed in this paper. The developed impedance is based on an asymmetric structured tree methodology as introduced and used for blood in Refs. ͓26-28͔. This cited work has been modified in order to represent the physiological environment in the pulmonary tree. The coupling of the 1D model to the 3D model is achieved utilizing a Dirichlet to Neumann approach introduced in Ref. ͓29͔ for arterial blood flow.

The aim of this paper is to apply this modeling approach to lung simulations and also to study the mechanics of the lung under artificial ventilation. In particular, using such an approach, the entire conducting airway can be modeled efficiently obtaining important pressure and flow information.

2 Methodology 2.1 3D Geometry Reconstruction. The geometry used in the simulations was segmented using the commercially available segmentation software MIMICS ͑Materialise, Leuven, Belgium͒. The standard CT scans had a resolution of 0.6 mm. This allowed up to seven generations of the bronchial tree to be obtained. Following segmentation the outlets of the geometry were cut normal to the flow direction and the geometry was exported in stereolithography ͑STL͒ format. Figure 1 shows the lung model with outlets cut normal to the flow direction.

2.2 Meshing. The geometry was then meshed in the commercially available meshing software HARPOON ͑Sharc, Manchester, UK͒. However, before meshing, the outlets here extruded within HARPOON in order to satisfy the requirements of the boundary conditions ͑see Sec. 2.4.2͒. In this software, a high quality mesh in the complex bronchial tree is a matter of minutes, an example mesh with extruded outlet extensions is shown in Fig. 2 . The final mesh contained approximately 2 ϫ 10 6 tetrahedral elements equating to 1.6ϫ 10 6 degrees of freedom.

The governing equations for airflow in the tracheobronchial region are the time dependent incompressible Navier-Stokes equations.

where u represents the velocity vector ͑u = ͓u , v , w͔͒, p is the air pressure, is the density, and is the dynamic viscosity. The numerical simulations were performed in our in-house multiphysics research code BACI. This code takes advantage of message passing interface ͑MPI͒ parallelization, hence, the solution of large scale simulations is achieved in a realistic time frame. The fluid domain was spatially discretised with stabilized finite elements. In time, the equations were discretised using a one-stepscheme. The resulting linear system was solved using the generalized minimal residual iterative solver with multilevel preconditioning, which is implemented in our code using the open-source package Aztec ͓30͔. For all simulations a time step size of 0.004s was used and the Newton iterations were considered converged when the residual dropped below 10 −6 .

The inflow boundary conditions used in the present simulations are based around normal breathing and mechanical ventilator prescribed flowrate under inspiratory conditions. For normal breathing, a sinusoidal inspiratory profile was utilized, q͑t͒ = A sin͑t͒, where A represents the flow amplitude and the circular frequency. This was done under two conditions with a mean flow of 0.2 l/s and 0.5 l/s, respectively, ͑Fig. 3͑a͒͒. This basic breathing model in the present study was only used as a benchmark since under nonmechanical ventilation conditions the upper airways ͑nasal pharynx͒ should also be included, in particular, because the glottis can produce turbulence in the trachea and the first few generations of the bronchial airways. Under mechanical ventilation conditions, the endotracheal tube is inserted into the level of the trachea, hence, we have not included the upper airways in our simulations.

For mechanical ventilation, a typical protective ventilator constant flow profile, see Fig. 3͑b͒ , was utilized. This profile is characterized by almost constant flow, apart from a rapid acceleration and deceleration period at the beginning and end of inspiration. Due to the acceleration periods and high frequency oscillations, the mechanical ventilation is expected to lead to different flow and pressure dynamics due to the considerably higher slew rates. The two aforementioned profiles ͑sinusoidal and constant flow venti-lation͒ cover the situations, which are available to the anaesthetist as control variables for ventilation.

At the walls of the model all velocity components are zeroed in order to enforce the no slip condition. Furthermore, since this study wants to focus on the effects of boundary conditions, the walls are assumed to be rigid. Future work will address the issue of FSI as has previously been investigated by our group ͓17͔, where simpler boundary conditions were used. 2.4.1 Impedance Tree. The airways beyond the 3D domain are modeled using a 1D approach previously described for arterial flow ͓26,27͔. This 1D approach is modified in order to model the pulmonary tree rather than the peripheral arterial tree. In order to calculate the impedance of the pulmonary tree, the airways are modeled as 1D flexible elastic tubes. To achieve this, the linearized 1D fluid momentum equation is used and is given by

where u is the velocity in the streamwise direction, t time, density, viscosity, and r the radius of the vessel. The approximation of linearized 1D momentum is reasonable for the pulmonary tree due to the high wave speed. In addition to the momentum equation, we also have the 1D continuity equation

where p is the pressure, q is the flowrate, and C is the vessel compliance, which is given by the following linearized state equation ͓22͔.

with r as the vessel radius, E as the elastic modulus, and h as the wall thickness. The wall thickness values are determined, based on vessel diameter, from the human Horsfield lung model given in Ref. ͓24͔ .

Equations ͑3͒ and ͑4͒ can then be formulated into an expression for impedance in the frequency domain by utilizing Womersley's solution, which is an analytical solution for oscillatory flow in a circular tube ͓31͔. The derived impedance equation expresses the impedance at the root of each airway segment ͑Z͑0,͒͒ as a function of the impedance in the downstream airway tree ͑Z͑L , ͒͒ and is given by

where g is the product of the wave speed c and the vessel compliance C, L is the vessel length, and the circular frequency. The wave speed c is given by

where A 0 is the cross-sectional area, F j is a function of the zeroth and first order Bessel functions as calculated from Womersley's solution, and is the density and compliance. The impedance of each individual branch is then summed in series and parallel in order to obtain the total impedance of a specific tree. At bifurcation points, the impedance of the parent vessel is related to the daughter vessels using a standard bifurcation condition.

where subscripts P, L, and R represent parent, left daughter, and right daughter vessels, respectively. In addition, scaling of the tree is achieved by introducing a relationship between parent and daughter vessels.

r L = ␣r P and r R = ␤r P ͑9͒

where ␣ and ␤ represent parameters for tree asymmetry. For example, if ␣ = ␤ = constant, the tree would represent a Weibel lung tree ͓15͔. In this paper, both symmetric and asymmetric trees were investigated. Asymmetric branching was based on averaged reported values in literature ͓32,33͔ with ␣ = 0.86 and ␤ = 0.7. Recursively applying Eqs. ͑6͒ and ͑8͒ from the root of the tree down to a critical radius provides the impedance of the lung tree in the frequency domain. In the current model, the impedance at the root of the tree is included but in all numerical examples in this study it was considered to be zero. Realistic nonzero impedance at the root of the tree will be addressed in future models, hence, taking into account peripheral compliance. This will be achieved by coupling to coupled multiscale simulations of the alveolar region ͓34,35͔. To calculate the time varying pressure at the root of a specific tree, the impedance must be first transformed into the time domain via

Finally, a time varying pressure is obtained at the root of the artificial tree by convoluting the impedance z͑t͒ with the history of flow at a specific outlet, i.e.,

Inflow curves utilized in the simulations "a… light breathing "Q mean =0.2 l/s… and light activity "Q mean =0.5 l/s… and "b… constant flow mechanical ventilation. Note the ventilation curve was chosen based on the high slew rates and oscillations, which occur in the profile to understand how these affect lung mechanics.

The above method was implemented into our in-house finite element computational fluid dynamics ͑CFD͒ code as part of a outflow boundary condition, which is explained in the following section.

Coupling. The coupling of the impedance model to the 3D domain was achieved using a previously described Dirichlet-Neumann approach, see Ref. ͓29͔ for full details. This method has the advantage of being mathematically well posed and has proven to be numerically stable. In addition, it also has the advantage that the solution of 1D and 3D domain is achieved simultaneously. The method was implemented both explicitly and implicitly into our finite element CFD code, however, little difference was observed between the methods and the explicit version was utilized due to the faster solution times. In brief, given the downstream domain pressure, p ‫ء‬ ͑t͒ ͑obtained from the impedance model͒, the time dependent pressure at the outflow of the 3D domain can be represented by

where ⌫ s is the outflow boundary surface, w is a weighting function, Ĩ is the identity tensor, is the viscous stress tensor, and n is the outward normal. The second term ͑͒ is assumed to be zero, which inherently assumes that the flow is fully developed ͑ =0͒.

To satisfy this condition, the outlets of the 3D model are extruded. All outlet extrusions were deemed suitably long such that the observed profile was parabolic and streamwise independent. Furthermore, the effect that these extensions have on the pressure is minimal. This was further verified using the traction free condition with and without the extensions and comparing the pressure at the equivalent location. Finally, the use of these pressure based boundary conditions, in particular, how the flow distributes, was verified by looking at the mean flow division ratio as a percentage of the inlet flow. The following distribution was obtained: right lobe ͑57%͒, left lobe ͑42%͒, right superior lobe ͑19%͒, right middle lobe ͑10%͒, right inferior lobe ͑28%͒, left superior lobe ͑24%͒, and left inferior lobe ͑18%͒. These flow divisions were deemed suitable based on reported literature values ͓16͔. The major advantage of this method is that no assumption needs to be made about the flow, hence, flow distribution is completely governed by the downstream region. Further validation was provided by considering the pressure drop calculated with this coupled approach, we considered the pressure drop over the whole model for the varying flowrates. This instantaneous drop was approximately 33 Pa for a flowrate of 0.5 l/s; if we consider the equivalent values reported in literature ͓36͔, the value shows good agreement. The bronchial pressure drop over the acceleratory phase is given in Fig. 4. 

The results of the simulations show well documented characteristics for flow in a bifurcating tree, in particular, centripetal acceleration due to curving vessels leads to flow phenomena such as flow recirculation and secondary flow. This all arises from the need to satisfy the continuity condition, leading to Dean like vortices forming in the daughter vessels. Furthermore, the nonuniform geometry causes acceleration of the fluid in the higher generations ͑for clarification the lowest generation is the trachea͒, which is not observed in more idealized geometries ͑see flow vectors in Fig. 5͒ .

The narrowing of the vessels in these locations leads to jetting of the flow and, hence, pressure decreases. The flow remains relatively disturbed throughout the model, reinforcing the need for 3D simulations in the lower generations of the bronchial tree. In the higher generations, the flow is relatively uniform meaning the application of the impedance boundary conditions in these loca-tions is acceptable.

The most notable change between standard traction free and impedance outflow conditions is the pressure distribution, both spatially and temporally. In Fig. 6 , the spatial pressure distribution is compared for the two different boundary conditions: traction free ͑up to seven generation model͒ and impedance ͑seven generation model plus up to 13 generations of the 1D tree depending on the outlet size͒. Evidently, there are marked differences. For the present model, the pressure variation between different outlets did not vary considerably due to the similarity in outlet size. The time history of pressure is completely different for the two different scenarios. From a physiological prospective the above is of utmost importance because determination and understanding of correct airway pressure drop plays an important role in lung mechanics. The above immediately reinforces the need to consider the downstream domain when modeling the tracheobronchial region.

To observe the difference in more detail, normalized values were plotted ͑Fig. 7͒. This shows that with the impedance conditions there is a slight narrowing of the pressure-flow loop; this indicates an increase in the resistive component of impedance.

Overall, the two loops are very similar exhibiting nonlinear behavior with peak pressure and flow basically coinciding. This is to be expected based on the similar outlet sizes throughout the model. The normalized pressure distribution within the 3D domain is effected only marginally by the impedance of the downstream domain. The result is rather interesting, as it suggests that the opposition to flow comes mainly from the upper parts of the airway rather than the impedance of the lower tree. However, the pressure change between impedance and traction free is significant and is likely to play a central role when fluid-structure interaction of the 3D domain is considered. Figure 8 shows the distribution of mean pressure along two different centerline paths: left lobe ͑gray͒ and right lobe ͑black͒. Overall there is a drop in the pressure descending the vessel. However, due to disturbed flow phenomena, pressure increases are observed in the higher generations. This result is in contrast to more idealized models, which lack the pressure variation due to local curvature changes. Hence the pressure in such models tends to drop in a relatively uniform manner along the vessel. The results here also indicate the strong dependence on asymmetry, which previously has been reported to be an important factor in the resulting pressure distribution. In particular the distribution between the left and right lung lobe bear no relationship. Furthermore, when the results are compared with the traction free model ͑data not shown͒ greater spatial gradients in pressure are observed. Finally, it is seen from the results that the largest pressure drop is observed in the higher generations.

The effect of different trees on the pressure distribution in the 3D domain is significant, however, this is to be expected as the impedance is dependent on the size of the tree, i.e., the smaller the tree the greater the impedance. This is demonstrated in Fig. 9 , where two different symmetric Weibel trees ͑␣ = ␤ = 0.5 and 0.9͒ are attached to the end of the model. Evidently the pressure in the 3D domain is considerable higher. In particular the resulting pressure in the faster diminishing tree is higher because essentially the tree can accommodate less flow ͑due to reduced cross-sectional area͒. This represents a constriction of the entire lung tree.

The above suggests that the depth of the tree is the most important factor, i.e., the number of generations. This means the structure of the downstream tree is relatively insignificant to the spatial pressure distribution, especially the difference between symmetric and asymmetric trees. Overall, the models utilizing an asymmetric tree are better because they are a closer representation of in vitro measured trees, however, some experimental validation based around pressure measurements is still required. In all presented models the effects on flow appear to be relatively minimal, there is a small change in the field, however, this appears to be relatively insignificant.

Finally, we consider how the boundary conditions can be used to simulate hypothetical disease in the human lung. The current impedance model has the advantage that no a priori assumptions need to be made about the flow since the flow in the 3D domain is governed by the peripheral vessels. This is a major advantage over previous impedance implementations. Here, we consider a light constriction of the left superior lobe, i.e., a lung of reduced functional residual capacity. This could represent clinical scenarios such as alectasis or a derecruited region of the lung. This model is only indicative and is not meant to reproduce the complex nature of true lung disease since lung disease formation tends to be heterogeneous. The constriction of the airways is modeled by increasing the impedance on the outlets. Here, the model responds by letting less flow into this lobe and naturally diverting flow to other regions of the lung. This shows that with our method we can simulate heterogeneous changes in the lung. Furthermore, consid- Fig. 6 Comparison of pressure contours, for light activity, at maximum inspiration: "a… traction free boundary condition "up to seven generations… and "b… impedance boundary condition"up to 17 generations…. With the impedance conditions the maximum pressure is approximately 44% higher. Note the lower limit of the scale is set to 5 Pa for visual comparison purposes and for the traction free condition the pressure at the outlets is virtually zero, meaning the pressure in the rest of the tree would be negative. ering Fig. 10 , the pressure-flow relationship is only really affected in the vicinity of the occlusion. In Figs. 10͑a͒ and 10͑b͒ , it is evident that there is a widening of the loop ͑hysteresis͒ indicating an increase in the phase difference between pressure and flow. However, in other parts of the lung the relationship remains very similar ͑between disease and undiseased scenarios͒. This suggests the response is quite local, i.e., pressure and flow is increased throughout the lung, however, the dynamic between the two remains similar.

Under mechanical ventilation complex flow dynamics result from both geometric and temporal influences. Ultimately geometry is the major factor in the resulting flow distribution, however, the rapid deceleration observed late in the inspiratory phase results in larger recirculation zones forming in the higher generations. The flowrate difference between the two different boundary conditions shows very similar profiles ͑data not shown͒. In addition, the instantaneous flowrate difference is maximal 4% higher under impedance conditions. The most no- Fig. 7 Normalized pressure-flow relationship in the bronchial airway. Evidently the pressure and flow are consistently in phase throughout the cycle. This is more noticeable for impedance conditions, where a narrowing of the loop is observed. For the impedance condition, the temporal development of pressure is more rapid. The difference in flowrate between the two conditions is marginal with the flowrate being overall reduced under impedance conditions "data not directly shown, however, at each point along the two above curves the flow is nearly identical…. table change between standard traction free and impedance outflow conditions is again the pressure distribution, both temporally and spatially. In particular the max pressure difference between the two different boundary conditions at peak inspiration is 27% in the second generation ͑traction free condition is 27% lower than impedance͒. The observed temporal gradients in pressure are interestingly high in the higher generations, however, the dynamics tended to be similar. This effect is driven by the nonuniform geometry.

The actual maximum pressure level for the first four to five generations generally varies only a small amount. Interestingly in the left lobe of the lung, the pressure over the cycle in the third generation is higher than in the second. This is driven by the flow route, more air is diverted into the left superior lobe, where the pressure is reduced ͑data not shown͒. This results in less air flowing into the left inferior lobe. The flow route is primarily driven by the curvature of the second generation of the left hand side of the lung, where the air is skewed toward the interior surface of the lung, hence flow continuity results in secondary currents, which lead to flow diverting back toward the anterior surface. This aforementioned phenomenon is primarily due to the nonuniformity of the geometry, meaning the cross-sectional area change plays an important role in the pressure and flow distribution. This particularly suggests that models based around circular cross sections fail to capture essential dynamics such as pressure variations due to cross-sectional changes. Figure 11 compares the pressure-flow relationship between ventilation and breathing. Clearly there is a change in this relationship. This is evidenced by the pressure-flow loop area increase, hence indicating that inertial effects are important ͑pressure is leading flow͒. This relationship is also seen to be much more nonlinear. This is specifically true during the acceleration phases. The small lag that does occur comes about due to transient fluid inertia, which results directly from higher slew rates observed during this ventilation protocol. Ultimately the flow in the lung can be characterized as quasi-steady ͑Womersley numbers ϳ1͒, which means the flow responds to the changing pressure gradient in a very rapid fashion. Both boundary conditions, under mechanical ventilation, exhibit phase differences; this indicates that the upper generations contribute more dominantly to the lung impedance than the peripheral vessels, this is also backed up by the dynamics being not significantly different. However, the peripheral vessels are obviously important for correct pressure distribution/level in the lower generations, as indicated by the large difference between pressure levels.

Again a hypothetical diseased lung is simulated. The results are shown in Fig. 12 . Most interesting is the result in the superior left lobe ͑Fig. 12͑b͒͒, where the increased impedance on the branch results in an in phase relationship ͑between pressure and flow͒ during acceleration and deceleration. This particularly shows a very different relationship to that of the healthy profile and highlights that the method of ventilation is important.

Modeling this reduced functional lung results in an overall increase in pressure in the lung, see Fig. 13 . In particular, the mean pressure plotted along a centerline ͑gray centerline from Fig. 8͒ is elevated throughout and at the outlet ͑located in the inferior left lobe͒ it is approximately 28% higher. The pressure is up to 38% higher in the vicinity of the occluded area. Similar observations are also seen in the right lobe. These results mean that the nondiseased part of the ventilated lung is potentially at risk from volutrauma due to an increase in lobular pressure over a normal functioning lung.

In this paper we have developed and used structured tree outflow boundary conditions for modeling the lower airways of the lung tree up to 17 generations. The coupling of this 1D model to the 3D domain uses a Dirichlet to Neumann approach. This approach is most general as it makes no assumption about the form of the downstream domain. Furthermore, the method is applied in a mathematically well posed formulation, i.e., avoiding a full Dirichlet problem. Currently the impedance model is based purely around a structured tree, however, the present implementation al- First the boundary conditions were compared under normal breathing conditions, demonstrating that there is a clear difference between the pressure levels. In addition, different trees were attached to the outlets, specifically, a Weibel tree with different scaling coefficients. This showed that the pressure level can change significantly. The faster scaling tree coefficient ͑␣ = ␤ = 0.5͒ represents a tree of reduced volume, causing elevated pressure in the lungs. This highlights the importance of the small airways on the pressure levels in the lower generations. Utilizing models in this way allows the effects of small airways to be investigated, which otherwise remains difficult ͓38͔.

We further investigated the effects of constant flow ventilation on lung mechanics. The high slew rate ventilation curve resulted in larger recirculation zones in higher generations, particularly during the deceleration phase. The higher frequencies of the mechanical ventilation profile led to a different pressure-flow dynamic, in particular during the rapid acceleration/deceleration phase and where the flow became slightly negative ͑compare Fig.  11͒ . The wider pressure-flow loop indicates a phase lag between pressure and flow; this is in contrast to normal breathing, where the pressure and flow remain essentially in phase. This different dynamic is due to the effects of transient fluid inertia. However, generally the flow in the lungs is quasi-steady since there is only a small lag between pressure and flow, i.e., flow responds rapidly to the changing pressure gradient. This is to be expected given the Womersley number is only marginally above one for the lower generations and below one for higher generations.

Overall, for all models the flowrate did not change considerably. The fact that the impedance does not change the flow magnitude significantly suggests that models, which are interested in flow related phenomena only, such as particle distribution, potentially do not require such complicated impedance type boundary conditions. At the very least, such simulations can consider only resistance at the outlet terminals. However, in making such an assumption, a significant number of generations must be considered in the model in order to capture the relevant pressure-flow dynamics, i.e., the impedance is predominantly from the upstream area of the model. This is supported by previous medical observations, where it is reported that around 10% of the resistance comes from the small airways ͓39,40͔. From a flow physics per-spective, this does not mean the number of generation levels, rather it requires that the Reynolds number at the outlet is sufficiently low.

Modeling diseased lungs is essential for understanding ventilation of diseased patients and can also provide insight into how different mechanical maneuvers affect the diseased lungs. Understanding how the lungs respond under mechanical ventilation is currently an unresolved question still trying to be elucidated by the medical community ͓41͔. Here, a partial obstruction of the superior left lobe was investigated in order to see the effects on the rest of the lung. This simulation also demonstrated the usefulness of the coupling method, where no assumption is made a priori about the flow into the peripheral airways, i.e., the flow is realistically governed by the impedance in this downstream region. The obstruction was not meant to model the complexity of lung disease but instead see the effects of the peripheral vessels on the lung dynamics. However, it could be indicative of a downstream derecruited region. The results demonstrated for both normal breathing and mechanical ventilation that the major changes in pressure-flow dynamics tended to be localized to the vicinity of the region, where the constricted vessels are located, i.e., the other regions of the lung showed increased pressure and flow, however, the relationship between the variables is virtually unchanged. In these other regions, the mean outlet pressure was on average approximately 28% higher, which could potentially result in overdistension ͑volutrauma͒ of the downstream alveoli due to elevated lobular pressure. The model here, although completely idealized, demonstrates that diseased lungs respond differently to healthy lungs. Further work is required in this area to understand heterogeneous lung disease. This model represents a step in the direction to understand the complex lung mechanics occurring during mechanical ventilation. One of the major areas, which requires understanding is the effect of different mechanical maneuvers on lung mechanics in a diseased state.

In literature, alternative approaches related to modeling the total tracheobronchial airways have been proposed ͓42,43͔. In Ref. ͓42͔ , a 15 generation model of the airway system was simulated by calculating the flow in a 3D Weibel triple bifurcation and then repeatedly reproducing this section in parallel, where the outflow from each section becomes the inflow for the next downstream bifurcation unit. In Ref. ͓43͔ , an idealized 17 generation model based on an anatomical data set was simulated under steady state conditions, however, their mesh used in the model was very course. Although these two methods provide some useful insights into flow structures in the conducting airways our approach has a number of advantages: We can model the whole conducting airway utilizing a realistic geometry for the first seven generations and then attach 1D trees as boundary conditions onto the 3D domain; the 1D domain takes into account the effects of wall movement; we are utilizing transient simulations; and we have a framework in which we have demonstrated that we can model diseased state conditions. Future work will also encompass the expiration phase, which is nontrivial from a computational perspective. The simulations will also take into account FSI effects using our vast experience in FSI and the efficient implementations in our in-house code ͓21͔. This is expected, given the change observed in pressure under impedance conditions to lead to different flow and stress dynamics in the tracheobronchial region. In particular significant changes are expected to be seen in the higher generations due to the reduced wall cartilage content. Currently the 1D model is attached to the outlet of the 3D model obtained from the CT data, future models will include artificial bifurcations based on morphological data to extend the 3D model a few extra generations, this is in order to limit any error associated with using a 1D model. Ultimately, we are working toward a whole lung model, hence, work is currently being undertaken to couple the model presented here to models of the alveolar region ͓44,34͔, therefore we will have a model of the tracheobronchial tree coupled to the acinar region. This is ex- pected to alter the results since this is a major part of the respiratory system. This is in the hope that further understanding of the lung, specifically the alveoli under healthy and diseased conditions can be ascertained. Finally, the present models have provided detailed insight into the effects of the peripheral airways on breathing mechanics and also the use of coupling for lung simulations.

Incidence and Outcomes of Acute Lung Injury

Effect of Mechanical Ventilation on Inflammatory Mediators in Patients With Acute Respiratory Distress Syndrome: A Randomized Controlled Trial

Mechanical Ventilation in Acute Lung Injury and Acute Respiratory Distress Syndrome

Effect of a Protective-Ventilation Strategy on Mortality in the Acute Respiratory Distress Syndrome

Low Mortality Associated With Volume Pressure Limited Ventilation With Permissive Hypercapnia in Severe Adult Respiratory Distress Syndrome

Higher Versus Lower Positive End-Expiratory Pressures in Patients With the Acute Respiratory Distress Syndrome

Mechanisms of Ventilator-Induced Lung Injury: The Clinician's Perspective

Acute Lung Injury

Protective Ventilation of Patients With Acute Respiratory Distress Syndrome

A Lung Computed Tomographic Assessment of Positive End-Expiratory Pressure-Induced Lung Overdistension

The Acute Respiratory Distress Syndrome

Modelling of Peak-Flow Wall Shear Stress in Major Airways of the Lung

Modeling the Bifurcating Flow in a Human Lung Airway

Airflow Structures and Nano-Particle Deposition in a Human Upper Airway Model

Morphometry of the Human Lung

Models of the Human Bronchial Tree

Fluid Structure Interaction in Lower Airways of CT-Based Lung Geometries

An Anatomically Based Hybrid Computational Model of the Human Lung and Its Application to Low Frequency Oscillatory Mechanics

Characteristics of the Turbulent Laryngeal Jet and Its Effect on Airflow in the Human Intra-Thoracic Airways

Characteristics of Airflow in a CT-Based Ovine Lung: A Numerical Study

Coupling Strategies for Biomedical Fluid-Structure Interaction Problems

Wave Propagation, Input Impedance, and Wall Mechanics of the Calf Trachea From 16 to 1,600 Hz

Airway Inhomogeneities Contribute to Apparent Lung Tissue Mechanics During Constriction

How Heterogeneous Bronchoconstriction Affects Ventilation Distribution in Human Lungs: A Morphometric Model

A Morphometric Model of Lung Mechanics for Time-Domain Analysis of Alveolar Pressures During Mechanical Ventilation

Structured Tree Outflow Condition for Blood Flow in Larger Systemic Arteries

Numerical Simulation and Experimental Validation of Blood Flow in Arteries With Structured-Tree Outflow Conditions

Fractal Network Model for Simulating Abdominal and Lower Extremity Blood Flow During Resting and Exercise Conditions

Outflow Boundary Conditions for Three-Dimensional Finite Element Modeling of Blood Flow and Pressure in Arteries

The Physics of Pulsatile Flow

CT-Based Geometry Analysis and Finite Element Models of the Human and Ovine Bronchial Tree

A Three-Dimensional Model of the Human Airway Tree

Towards a Comprehensive Computational Model for the Respiratory System

A Nested Dynamic Multi-Scale Approach for 3D Problems Accounting for Micro-Scale Multi-Physics

Pulmonary Fluid Dynamics

Generation of an Anatomically Based Three-Dimensional Model of the Conducting Airways

The Role of Small Airways in Lung Disease

Small Airways: A Time To Revisit?

The Physiology of Small Airways

Alveolar Recruitment in Acute Lung Injury

An Adjustable Triple-Bifurcation Unit Model for Air-Particle Flow Simulations in Human Tracheobronchial Airways

Computational Model of Airflow in Upper 17 Generations of Human Respiratory Tract

Modeling the Mechanical Behavior of Lung Tissue at the Micro-Level

We acknowledge the support of Deutsche Forschungsgemeinschaft ͑DFG͒ through Project No. WA-1521/8. The CT images of the tracheobronchial tree have been provided by Dr. med Ulrich Barkow and Dr. med. univ. Andreas Strohmaier from "Diagnostische Radiologie," Stuttgart, Germany.