- Email: [email protected]
Respiratory flow in obstructed airways
ARTICLE IN PRESS
Journal of Biomechanics 39 (2006) 2743–2751 www.elsevier.com/locate/jbiomech www.JBiomech.com
Respiratory flow in obstructed airways X.L. Yang, Yang Liu, H.Y. Luo Department of Mechanical Engineering, The Hong Kong Polytechnic University, Hung Hom, Kowloon, Hong Kong Accepted 10 October 2005
Abstract Chronic obstructive pulmonary disease (COPD) is one of the most common diseases in human community. The COPD always results in inflammation that leads to narrowing and obstruction of the airways. The obstructive airways have significant effect on respiratory flow. In order to understand the flow phenomenon in such obstructive airways, four three-dimensional four-generation lung models based on the 23-generation model of Weibel [1963. Morphometry of the Human Lung. Springer, Academic Press, Berlin, New York] are generated. The fully three-dimensional incompressible laminar Navier–Stokes equations are solved using computational fluid dynamics (CFD) solver on unstructured tetrahedral meshes. Therein, a symmetric four-generation airway model is served as the reference, the other three models are considered to be obstructed at each generation, respectively. The calculation results show that the obstructive airway has significant influence on the air flow in both up- and down-stream airways and it even results in flow separation in the conjunction region. The re-circulation cell blocks the air from entering the downstream branches. This may be the reason why COPD patients should breathe gently, and this also provides some valuable information for medicine powder deposition. r 2005 Elsevier Ltd. All rights reserved. Keywords: Respiratory flow; COPD; Four-generation; CFD
1. Introduction Chronic obstructive pulmonary disease (COPD) is one of the most common diseases in the world and it is caused by blocking of the airways in the lung which invariably results from heavy smoking and inhalation of pollutants of various kinds. Numerous studies have been carried out on this topic and the typical model is a symmetric bifurcation airway. However, chronic bronchitis is an inflammation of the bronchi which alters the branching configuration significantly. The inflamed bronchi are no longer symmetric or regular asymmetric airways. From previous studies (Liu et al., 2002, 2003), the flow rate and secondary flow pattern in human lung are quite sensitive to the tube diameter and bifurcation configuration. The flow characteristics of the asymmetric airway are totally Corresponding author. Tel.:+852 2766 7814; fax: +852 2365 4703.
E-mail address: [email protected] (Y. Liu). 0021-9290/$ - see front matter r 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.jbiomech.2005.10.009
different from that of the symmetric one. Therefore, obstructions in airways could significantly alter the velocity distribution, air flow rate, and pressure distribution. Most of the early studies focused on the flow and deposition in a symmetric bifurcation model under symmetric flow conditions. From these studies (Schroter and Sudlow, 1969; Pedley et al., 1977; Isabey and Chang, 1981; Snyder and Olson, 1989; Zhao and Lieber, 1994), skewed velocity profiles and two symmetric eddies in the two-generation experimental bifurcation flow models were found. Detailed air/particle transport phenomena as well as local deposition patterns and surface densities of deposited particle in bifurcating airways are difficult to obtain experimentally. Consequently, computational fluid dynamics (CFD) simulation is an alternative and an effective tool to acquire such information and enhance understanding of this difficult problem. Several multi-dimensional simulations have been carried out to
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
2744
G8-4
generation model of Weibel (1963). The diameter of each generation is equal to that of the fifth to eighth in the model of Weibel (1963), respectively. The bifurcation angle is 701. Detailed geometric parameters are tabulated in Table 1. In order to study the effect of COPD on bifurcation flow, three additional models are individually modified from the symmetric model (model 1). For model 2, the diameter of one of the second generation tubes is smoothly contracted to half of the original diameter (i.e., the diameter of the model 1, D6) at the middle of this airway. For model 3, the diameter of one of the lateral third generation tubes is reduced to the half of its original size; and for model 4, the diameter of one of the medial third generation tubes is reduced to the half of its original size. The fully three-dimensional incompressible laminar Navier–Stokes equations are solved using CFD solver based on unstructured meshes (Liu et al., 2002). The governing equations are solved sequentially using the segregated method. The convective terms are all discretized using second order upwind scheme and the SIMPLEC method is used for the pressure–velocity coupling. In addition, the pressure interpolation scheme is specified as second order. The multigrid method is invoked to increase the speed of convergence of numerical solutions. To ensure the solution accurate and stable, conservative under-relaxation factors are selected to be 0.3 for pressure and 0.5 for momentum. A parabolic velocity profile is imposed on the inlet; the static pressure is set to be zero at all outlets; no-slip wall boundary condition is imposed on all solid walls. Five Reynolds numbers, which are based on the mean velocity on the inlet and the diameter of the first generation tube, ranging from 300 to 1500 with increment of 300, are considered in each model. The cell numbers for the four computational models are 948285, 956559, 974047 and 975388, respectively. These numbers are determined by using different meshes, from coarse to progressively fine, until the calculated mass flow rates are mesh-convergent to within a prescribed tolerance (0.5%). A refined mesh has been employed near the walls and the conjunctions where the velocity gradient may be larger. Details of the algorithm have been discussed elsewhere, the computer simulation model for airflow has
G8-3
Table 1 The geometric parameters of the models
study bifurcation flow and particle transport in the upper respiratory tract, the bronchial airways and the acinar region of the lung (Martonen et al., 1994; Balashazy et al., 1999; Lee et al., 2000; Comer et al., 2001; Darquenne, 2001; Liu et al., 2002; Hegedus et al., 2004; and Miguel et al., 2004). These studies have provided additional detailed information of particle behavior in the lung airways and showed that this behavior is inherently linked to the fluid flow patterns within the airways. However, all the studies only considered a straight smooth airway, i.e. the effect of obstruction was ignored. In this study, we are concerned with the threedimensional flow in four-generation obstructed bifurcation airways. A symmetric four-generation airway model is served as the reference, and the other three models are considered to be obstructed at either the second generation or the third generation airways. The objective is to investigate the influence of obstructed airway on the flow pattern, air flow rate and pressure drop.
2. Numerical method A schematic view of the four-generation models adopted in this study is shown in Fig. 1. The bold solid lines indicate a symmetric model, which is constructed based on the fifth to eighth generations of the 23-
G8-8 G8-7 x G7-4 z
G8-6
y
G7-3 G6-2
L5 L8
G8-5
G6 -1 G5
L7 L6
G7-2 G7-1
G8-2 G8-1 Fig. 1. The schematic view of the computational model.
Generation number
Diameter D (mm)
Length L (mm)
5 6 7 8
3.5 2.8 2.3 1.86
10.7 9.0 7.6 6.4
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
2745
Fig. 2. The flow patterns of model 1: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) path lines of Re ¼ 1500.
been validated with experimental velocity data of Zhao and Lieber (1994) for a two-generation symmetric bifurcation (Liu et al., 2002); agreement between experimental and numerical results is excellent.
3. Results and discussion To thoroughly investigate the effect of COPD on bifurcation flow, five Reynolds numbers, 300, 600, 900, 1200 and 1500, are selected to simulate the respiratory flow, where the Reynolds number is defined as Re ¼ rUD=m, in which r is the air density, U is the mean inlet velocity, D is the diameter of the inlet tube, and m is the dynamic viscosity.
3.1. Flow pattern To discuss the effect of COPD on the flow pattern comparatively, the symmetric model (model 1) is considered first. Fig. 2 shows the flow vectors of the model 1 (symmetric model) in the bifurcation plane for Re ¼ 300 and 1500, respectively. Airflow splits into branches G6-1 and G6-2 at the carina ridge. Centrifugal force pushes more air to move toward inner walls of branches G6-1 and G6-2, consequently a skewed velocity profile is clearly presented. The skewed velocity distribution makes more air flows into the two medial branches of the third generation, G7-2 and G7-3. With increasing Re, the skewness becomes significant and more air flows into G7-2 and G7-3. To see clearly the flow rate distribution, the path lines at Re ¼ 1500 are
ARTICLE IN PRESS 2746
X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
plotted in Fig. 2(c). The flow rate distribution is exactly symmetric about G6-1 and G6-2, and there is no recirculation or stagnation throughout the airways. The airflow strongly depends on the airway configuration, when one of the second generation airways is obstructed (model 2), it inevitably alters the flow pattern. Fig. 3 shows the flow patterns at different Re for model 2 in which G6-1 is obstructed. Generally, the contracted tube has altered the flow field significantly, particularly for the downstream airways. The contracted airway (G6-1) behaves like a ‘‘jet’’, this jet increases the flow momentum, and the skewed velocity profile does not have enough length to recover before the second conjunction, consequently most of the fluid enters the medial branch (G7-2) and the flow vectors in the lateral branch (G7-1) is very weak. With increasing Re (increasing momentum), the ‘‘jet’’ effect is even stronger,
and the flow is separated after the throat. Under the effect of the centrifugal force, most of the high momentum fluid flows into the medial branch (G7-2), and the re-circulation occurs in the lateral branch (G71), as shown in the enlarged view in Fig. 3(c). This recirculation cell blocks the fluid from entering the lateral branch. This may explain why the deep breath usually does not help for the COPD patients. The COPD patients should breathe gently and slowly, then the oxygen can reach every alveolus in the lung. If they take a deep breath, the re-circulation downstream of the obstructed part may block the oxygen from reaching the downstream alveoli. To show the airflow distribution thoroughly, the path line is also calculated and indicated in the bifurcation plane for Re ¼ 1500 in Fig. 3(d). Generally, due to the obstructive effect in branch G6-1, most of the fluid flows into branch G6-2, and the medial
Fig. 3. The flow patterns of model 2: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) enlarged view; (d) path lines of Re ¼ 1500.
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
branch G7-3 takes more fluid because of the centrifugal effect. In the obstructed branch G6-1, most of the fluid goes into the medial branch G7-2. The flow skews toward the inner wall and a much smaller portion of the fluid enters the lateral branch G7-1. The jet flow in the obstruction generates a strong re-circulation at the entrance of the lateral branch G7-1 and a weak recirculation in the mid-airway of G7-2. The re-circulation in G7-1 induces another much weaker re-circulation at the entrance of G8-2. Due to the obstruction in the midairway of G6-1, a slight separation occurs at the entrance of G6-1, consequently a weak stagnation area forms there. The obstruction in G6-1 generates several re-circulations in daughter-branches which block seriously the fluid from entering downstream airways. In current four-generation configuration, the obstruction in G6-1 can help us to understand the effect of
2747
obstruction on downstream airways, but it is difficult to evaluate the influence on upstream bifurcation airways. To evaluate the influence of obstruction on both up- and down-stream airways, both models 3 and 4 are obstructed in the third generation airway. Fig. 4 shows the flow patterns of model 3 in which the lateral airway G7-1 is obstructed. Due to the resistance of the obstructed airway in G7-1 and the effect of the centrifugal force in G6-1, most of the fluid in G6-1 flows into the medial branch G7-2. With increasing Re, the velocity profile in G6-1 skews significantly towards to the medial branch G7-2. Since the flow is weak in lateral branch G7-1, the contracted tube makes the velocity profile almost symmetric after the ‘‘throat’’ which can be seen from the enlarged view in Fig. 4(c). The recirculation exists at both the entrance of G8-1 and G8-2, even the re-circulation is very weak. The path lines in
Fig. 4. The flow patterns of model 3: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) enlarged view; (d) path lines of Re ¼ 1500.
ARTICLE IN PRESS 2748
X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
Fig. 4(d) clearly show that the obstruction in G7-1 not only generates a pair of almost symmetric re-circulation in daughter airways, but also results in a stagnation region in upstream G6-1. This indicates that, in bifurcating airways, the obstruction affects the flow in both down- and up-stream airways. Fig. 5 shows the flow patterns in the bifurcation plane for model 4 in which the G7-2 airway is contracted. With increasing Re, the flow vectors are skewed significantly towards the inner wall. The obstructed airway functions like a ‘‘jet’’, for low momentum jet flow (Re ¼ 300), the velocity profile is almost symmetric at the outlets G8-3 and G8-4; for high momentum jet flow (Re ¼ 1500), the velocity profile is skewed towards the inner wall, indicating that the flow does not have enough length to recover due to the high momentum. Because of the resistance of the obstructed airway, most of the fluid
in G6-1 flows into the lateral branches G7-1, even the velocity profile skews towards the medial branch G7-2 in the second generation branch G6-1. From the path line in Fig. 5(d), most of the jet fluid flows into G8-3, and the obstruction generates the re-circulation both down- and up-stream airways. Particularly, the stagnation in up-stream G6-1 locates at almost the same position as that of model 3. 3.2. Mass flow partition In a health lung airway, the air should spread uniformly to each alveolus. However, in an obstructed airway, the air mass flow rate is biased at the first conjunction due to the non-uniform resistance downstream. The imbalance of the mass flow affects the respiration extremely. Therefore, the mass flow partition
Fig. 5. The flow patterns of model 4: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) enlarged view; (d) path lines of Re ¼ 1500.
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
is crucial to understand the effect of COPD on the respiratory flow. In this study, the ratio of mass flow rate of branches i to j is defined as _i m , _j m
0.8
(1)
where m ¯ i;j is the ratio of mass flow rate of branches i to j, i and j indicate the branch number, e.g. 5, 6-1, 6-2, 7-1 _ i and m _ j are the mass flow rates through branches etc., m i and j, respectively. Fig. 6 shows the influence of the obstruction on the flow partition in the obstructed branches, i.e. the flow rate in G6-1 to the total flow rate in G5. In symmetric airway (model 1), the air flow rate in G6-1 takes half of the total flow rate, indicating a balanced or symmetric flow rate between the two big branches, G6-1 and G6-2. When the second generation G6-1 is obstructed, the flow rate in G6-1 reduces significantly, and the ratio is only about 25%. The obstruction in the third generation (G71 and G7-2) results in a reduction of flow rate in branch G6-1, and the influence of the lateral branch (G7-1) and medial branch (G7-2) is almost identical. For the air flow ratio between the third generation (G7-1) and second generation (G6-1), as shown in Fig. 7, the Re may have significant influence on the flow partition. In the symmetric airway (model 1), the lateral branch G7-1 takes slightly less air compared to the medial branch G7-2 due to the centrifugal effect. When G6-1 is obstructed (model 2), the air flow ratio decreases significantly with increasing Re. It is not surprising, because the high momentum fluid induces stronger recirculation at entrance of G7-1, which blocks the fluid from entering it. It is understandable that, a small amount of air flows into G7-1 when the third generation
0.5
m6-1,5
0.4 Model 1 Model 2 Model 3 0.3
Model 1 Model 2 Model 3 Model 4
0.6
m7-1,6 -1
m ¯ i;j ¼
2749
0.4
0.2
500
1000
1500
Re Fig. 7. The variation of mass flow ratio of m ¯ 71;61 with Re.
G7-1 is obstructed; and most of the air enters G7-1 when G7-2 is obstructed. Air flow ratio is the effective way to evaluate the influence of fluid momentum on the re-circulation blockage. Figs. 8(a) and (b) show the variation of m ¯ 81;71 and m ¯ 84;72 with Re, respectively. Both the G8-1 and G8-4 airways are lateral branches for G6-1 branch, if the velocity profile in G6-1 is symmetric, i.e., it is either uniform or parabolic, the m ¯ 81;71 and m ¯ 84;72 should be the same in the symmetric airways. However, due to the centrifugal effect in G6-1 and obstruction in airways, the m ¯ 81;71 and m ¯ 84;72 are totally different and strongly dependent on Re. When G6-1 is obstructed (model 2), with increasing Re, m ¯ 81;71 increases greatly, indicating that the re-circulation at the entrance of G8-2 becomes stronger (Refer to Fig. 3(d)) which blocks the air from entering G8-2; with increasing momentum, the m ¯ 84;72 decreases significantly, the reason is the re-circulation in the mid-airway of G7-2 becomes stronger and blocks the air from getting into the G8-4 airway (Refer to Fig. 3(d)). When G7-1 is obstructed (model 3), both m ¯ 81;71 and m ¯ 84;72 increase with increasing Re. When G7-2 is obstructed (model 4), the momentum has little influence on m ¯ 81;71 , but m ¯ 84;72 drops dramatically with increasing Re.
Model 4
3.3. Pressure drop behavior
0.2
500
1000
1500
Re
Fig. 6. The variation of mass flow ratio of m ¯ 61;5 with Re.
The pressure drop in lung airways is the driving force of the respiratory flow, in the present calculation, the inlet mass flow rate is fixed and the outlet pressure is constant, therefore the obstructed airway may result in different flow resistance. The pressure drop coefficient
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
2750
Model 1 Model 2 Model 3 Model 4
Model 1 Model 2 Model 3 Model 4
3
Cp
m8 -1,7 -1
0.6
2 0.5
1 500
1000
(a)
1500
500
Re
1000
1500
Re Fig. 9. The variation of pressure drop coefficient C P with Re.
the biggest and this results in a much higher flow resistance, particularly for the low momentum fluid. This also explains why the COPD patients feel difficult to take breath.
m8 -4,7 -2
0.5
4. Conclusions 0.4
The effect of COPD on respiratory flow is numerically studied on four different four-generation models. The fully three-dimensional incompressible laminar Navier– Stokes equations are solved by a CFD solver. The numerical calculation leads to the following conclusions:
Model 1 Model 2 Model 3 Model 4 0.3 500
(b)
1000
1500
Re
Fig. 8. The variation of mass flow ratio of m ¯ 81;71 and m ¯ 84;72 with Re: (a) m ¯ 81;71 ; (b) m ¯ 84;72 .
C P is defined as CP ¼
Dptotal pin;dynamic
,
(2)
where Dptotal is the mass-weighted average total pressure drop from inlet to outlets, pin;dynamic is mass-weighted average dynamic pressure at the inlet. The variation of C P with Re for different model is plotted in Fig. 9. Generally, the pressure drop coefficient C P decreases greatly with increasing Re, it is not surprising since it is consistent with the Bernoulli’s equation. The C P distribution for models 1, 3 and 4 are quite similar, but the C P of model 2 is much higher than that of other models. The most likely explanation of this phenomenon is that the diameter reduction of model 2 is
(1) The obstructed airway alters the flow field significantly, a strong separation region exists behind the ‘‘throat’’ at higher Reynolds number. In a bifurcation airway, the obstruction may generate recirculation both upstream and downstream. These re-circulation cells block the air from entering the downstream branches which indicates that a deep breath may block the oxygen from reaching the alveoli for a COPD patient. (2) The obstructed airway has significant influence on the daughter branches and it increases the flow resistance significantly.
Acknowledgements Support given by the Research Grants Council of the Government of the HKSAR under Grant No. PolyU 5273/04E PolyU and by the Hong Kong Polytechnic University under Central Research Grant Nos. A-PG08 and G-T677 is gratefully acknowledged.
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
References Balashazy, I., Hofmann, W., Heistracher, T., 1999. Computation of local enhancement factors for the quantification of particle deposition patterns in airway bifurcations. Journal of Aerosol Science 30, 185–203. Comer, J.K., Kleinstreuer, C., Zhang, Z., 2001. Flow structures and particle deposition patterns in double-bifurcation airway models, Part 1, air flow fileds. Journal of Fluid Mechanics 435, 25–54. Darquenne, C., 2001. Aerosol transport in the lung. In: Prisk, G.K., Paiva, M., West, J.B. (Eds.), Gravity and the Lung: Lessons from Microgravity. Marcel Dekker, Inc., New York, pp. 117–148. Hegedus, C.J., Balashazy, I., Farkas, A., 2004. Detailed mathematical description of the geometry of airway bifurcations. Respiratory Physiology & Neurobiology 141, 99–114. Isabey, D., Chang, H.K., 1981. Steady and unsteady pressure–flow relationships in central airways. Journal of Applied Physiology 51, 1338–1348. Lee, J.W., Lee, D.Y., Kim, W.S., 2000. Dispersion an aerosol bolus in a double bifurcation. Journal of Aerosol Science 31, 491–505. Liu, Y., So, R.M.C., Zhang, C.H., 2002. Modeling the bifurcating flow in a human lung airway. Journal of Biomechanics 35, 465–473.
2751
Liu, Y., So, R.M.C., Zhang, C.H., 2003. Modeling the bifurcating flow in an asymmetric human lung airway. Journal of Biomechanics 35, 951–959. Martonen, T.B., Guan, X., Schreck, R.M., 1994. Effects of carinal ridge shapes on lung airstreams. Aerosol Science and Technology 21, 119–136. Miguel, A.F., Reis, A.H., Aydin, M., 2004. Aerosol particle deposition and distribution in bifurcating ventilation ducts. Journal of Hazardous Materials B 116, 249–255. Pedley, T.J., Schroter, R.C., Sudlow, M.F., 1977. Gas flow and mixing in the airways. In: West, J. (Ed.), Bioengineering Aspects of the Lung. Marcel Dekker, NewYork, pp. 163–265. Schroter, R.C., Sudlow, M.F., 1969. Flow patterns in models of the human bronchial airways. Respiratory Physiology 7, 341–355. Snyder, B., Olson, D.E., 1989. Flow development in a model airway bronchus. Journal of Fluid Mechanics 207, 378–392. Weibel, E.R., 1963. Morphometry of the Human Lung. Springer, Academic Press, Berlin, New York. Zhao, Y., Lieber, B.B., 1994. Steady inspiratory flow in a model symmetric bifurcation. ASME Journal of Biomechanical Engineering 116, 488–496.
Journal of Biomechanics 39 (2006) 2743–2751 www.elsevier.com/locate/jbiomech www.JBiomech.com
Respiratory flow in obstructed airways X.L. Yang, Yang Liu, H.Y. Luo Department of Mechanical Engineering, The Hong Kong Polytechnic University, Hung Hom, Kowloon, Hong Kong Accepted 10 October 2005
Abstract Chronic obstructive pulmonary disease (COPD) is one of the most common diseases in human community. The COPD always results in inflammation that leads to narrowing and obstruction of the airways. The obstructive airways have significant effect on respiratory flow. In order to understand the flow phenomenon in such obstructive airways, four three-dimensional four-generation lung models based on the 23-generation model of Weibel [1963. Morphometry of the Human Lung. Springer, Academic Press, Berlin, New York] are generated. The fully three-dimensional incompressible laminar Navier–Stokes equations are solved using computational fluid dynamics (CFD) solver on unstructured tetrahedral meshes. Therein, a symmetric four-generation airway model is served as the reference, the other three models are considered to be obstructed at each generation, respectively. The calculation results show that the obstructive airway has significant influence on the air flow in both up- and down-stream airways and it even results in flow separation in the conjunction region. The re-circulation cell blocks the air from entering the downstream branches. This may be the reason why COPD patients should breathe gently, and this also provides some valuable information for medicine powder deposition. r 2005 Elsevier Ltd. All rights reserved. Keywords: Respiratory flow; COPD; Four-generation; CFD
1. Introduction Chronic obstructive pulmonary disease (COPD) is one of the most common diseases in the world and it is caused by blocking of the airways in the lung which invariably results from heavy smoking and inhalation of pollutants of various kinds. Numerous studies have been carried out on this topic and the typical model is a symmetric bifurcation airway. However, chronic bronchitis is an inflammation of the bronchi which alters the branching configuration significantly. The inflamed bronchi are no longer symmetric or regular asymmetric airways. From previous studies (Liu et al., 2002, 2003), the flow rate and secondary flow pattern in human lung are quite sensitive to the tube diameter and bifurcation configuration. The flow characteristics of the asymmetric airway are totally Corresponding author. Tel.:+852 2766 7814; fax: +852 2365 4703.
E-mail address: [email protected] (Y. Liu). 0021-9290/$ - see front matter r 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.jbiomech.2005.10.009
different from that of the symmetric one. Therefore, obstructions in airways could significantly alter the velocity distribution, air flow rate, and pressure distribution. Most of the early studies focused on the flow and deposition in a symmetric bifurcation model under symmetric flow conditions. From these studies (Schroter and Sudlow, 1969; Pedley et al., 1977; Isabey and Chang, 1981; Snyder and Olson, 1989; Zhao and Lieber, 1994), skewed velocity profiles and two symmetric eddies in the two-generation experimental bifurcation flow models were found. Detailed air/particle transport phenomena as well as local deposition patterns and surface densities of deposited particle in bifurcating airways are difficult to obtain experimentally. Consequently, computational fluid dynamics (CFD) simulation is an alternative and an effective tool to acquire such information and enhance understanding of this difficult problem. Several multi-dimensional simulations have been carried out to
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
2744
G8-4
generation model of Weibel (1963). The diameter of each generation is equal to that of the fifth to eighth in the model of Weibel (1963), respectively. The bifurcation angle is 701. Detailed geometric parameters are tabulated in Table 1. In order to study the effect of COPD on bifurcation flow, three additional models are individually modified from the symmetric model (model 1). For model 2, the diameter of one of the second generation tubes is smoothly contracted to half of the original diameter (i.e., the diameter of the model 1, D6) at the middle of this airway. For model 3, the diameter of one of the lateral third generation tubes is reduced to the half of its original size; and for model 4, the diameter of one of the medial third generation tubes is reduced to the half of its original size. The fully three-dimensional incompressible laminar Navier–Stokes equations are solved using CFD solver based on unstructured meshes (Liu et al., 2002). The governing equations are solved sequentially using the segregated method. The convective terms are all discretized using second order upwind scheme and the SIMPLEC method is used for the pressure–velocity coupling. In addition, the pressure interpolation scheme is specified as second order. The multigrid method is invoked to increase the speed of convergence of numerical solutions. To ensure the solution accurate and stable, conservative under-relaxation factors are selected to be 0.3 for pressure and 0.5 for momentum. A parabolic velocity profile is imposed on the inlet; the static pressure is set to be zero at all outlets; no-slip wall boundary condition is imposed on all solid walls. Five Reynolds numbers, which are based on the mean velocity on the inlet and the diameter of the first generation tube, ranging from 300 to 1500 with increment of 300, are considered in each model. The cell numbers for the four computational models are 948285, 956559, 974047 and 975388, respectively. These numbers are determined by using different meshes, from coarse to progressively fine, until the calculated mass flow rates are mesh-convergent to within a prescribed tolerance (0.5%). A refined mesh has been employed near the walls and the conjunctions where the velocity gradient may be larger. Details of the algorithm have been discussed elsewhere, the computer simulation model for airflow has
G8-3
Table 1 The geometric parameters of the models
study bifurcation flow and particle transport in the upper respiratory tract, the bronchial airways and the acinar region of the lung (Martonen et al., 1994; Balashazy et al., 1999; Lee et al., 2000; Comer et al., 2001; Darquenne, 2001; Liu et al., 2002; Hegedus et al., 2004; and Miguel et al., 2004). These studies have provided additional detailed information of particle behavior in the lung airways and showed that this behavior is inherently linked to the fluid flow patterns within the airways. However, all the studies only considered a straight smooth airway, i.e. the effect of obstruction was ignored. In this study, we are concerned with the threedimensional flow in four-generation obstructed bifurcation airways. A symmetric four-generation airway model is served as the reference, and the other three models are considered to be obstructed at either the second generation or the third generation airways. The objective is to investigate the influence of obstructed airway on the flow pattern, air flow rate and pressure drop.
2. Numerical method A schematic view of the four-generation models adopted in this study is shown in Fig. 1. The bold solid lines indicate a symmetric model, which is constructed based on the fifth to eighth generations of the 23-
G8-8 G8-7 x G7-4 z
G8-6
y
G7-3 G6-2
L5 L8
G8-5
G6 -1 G5
L7 L6
G7-2 G7-1
G8-2 G8-1 Fig. 1. The schematic view of the computational model.
Generation number
Diameter D (mm)
Length L (mm)
5 6 7 8
3.5 2.8 2.3 1.86
10.7 9.0 7.6 6.4
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
2745
Fig. 2. The flow patterns of model 1: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) path lines of Re ¼ 1500.
been validated with experimental velocity data of Zhao and Lieber (1994) for a two-generation symmetric bifurcation (Liu et al., 2002); agreement between experimental and numerical results is excellent.
3. Results and discussion To thoroughly investigate the effect of COPD on bifurcation flow, five Reynolds numbers, 300, 600, 900, 1200 and 1500, are selected to simulate the respiratory flow, where the Reynolds number is defined as Re ¼ rUD=m, in which r is the air density, U is the mean inlet velocity, D is the diameter of the inlet tube, and m is the dynamic viscosity.
3.1. Flow pattern To discuss the effect of COPD on the flow pattern comparatively, the symmetric model (model 1) is considered first. Fig. 2 shows the flow vectors of the model 1 (symmetric model) in the bifurcation plane for Re ¼ 300 and 1500, respectively. Airflow splits into branches G6-1 and G6-2 at the carina ridge. Centrifugal force pushes more air to move toward inner walls of branches G6-1 and G6-2, consequently a skewed velocity profile is clearly presented. The skewed velocity distribution makes more air flows into the two medial branches of the third generation, G7-2 and G7-3. With increasing Re, the skewness becomes significant and more air flows into G7-2 and G7-3. To see clearly the flow rate distribution, the path lines at Re ¼ 1500 are
ARTICLE IN PRESS 2746
X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
plotted in Fig. 2(c). The flow rate distribution is exactly symmetric about G6-1 and G6-2, and there is no recirculation or stagnation throughout the airways. The airflow strongly depends on the airway configuration, when one of the second generation airways is obstructed (model 2), it inevitably alters the flow pattern. Fig. 3 shows the flow patterns at different Re for model 2 in which G6-1 is obstructed. Generally, the contracted tube has altered the flow field significantly, particularly for the downstream airways. The contracted airway (G6-1) behaves like a ‘‘jet’’, this jet increases the flow momentum, and the skewed velocity profile does not have enough length to recover before the second conjunction, consequently most of the fluid enters the medial branch (G7-2) and the flow vectors in the lateral branch (G7-1) is very weak. With increasing Re (increasing momentum), the ‘‘jet’’ effect is even stronger,
and the flow is separated after the throat. Under the effect of the centrifugal force, most of the high momentum fluid flows into the medial branch (G7-2), and the re-circulation occurs in the lateral branch (G71), as shown in the enlarged view in Fig. 3(c). This recirculation cell blocks the fluid from entering the lateral branch. This may explain why the deep breath usually does not help for the COPD patients. The COPD patients should breathe gently and slowly, then the oxygen can reach every alveolus in the lung. If they take a deep breath, the re-circulation downstream of the obstructed part may block the oxygen from reaching the downstream alveoli. To show the airflow distribution thoroughly, the path line is also calculated and indicated in the bifurcation plane for Re ¼ 1500 in Fig. 3(d). Generally, due to the obstructive effect in branch G6-1, most of the fluid flows into branch G6-2, and the medial
Fig. 3. The flow patterns of model 2: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) enlarged view; (d) path lines of Re ¼ 1500.
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
branch G7-3 takes more fluid because of the centrifugal effect. In the obstructed branch G6-1, most of the fluid goes into the medial branch G7-2. The flow skews toward the inner wall and a much smaller portion of the fluid enters the lateral branch G7-1. The jet flow in the obstruction generates a strong re-circulation at the entrance of the lateral branch G7-1 and a weak recirculation in the mid-airway of G7-2. The re-circulation in G7-1 induces another much weaker re-circulation at the entrance of G8-2. Due to the obstruction in the midairway of G6-1, a slight separation occurs at the entrance of G6-1, consequently a weak stagnation area forms there. The obstruction in G6-1 generates several re-circulations in daughter-branches which block seriously the fluid from entering downstream airways. In current four-generation configuration, the obstruction in G6-1 can help us to understand the effect of
2747
obstruction on downstream airways, but it is difficult to evaluate the influence on upstream bifurcation airways. To evaluate the influence of obstruction on both up- and down-stream airways, both models 3 and 4 are obstructed in the third generation airway. Fig. 4 shows the flow patterns of model 3 in which the lateral airway G7-1 is obstructed. Due to the resistance of the obstructed airway in G7-1 and the effect of the centrifugal force in G6-1, most of the fluid in G6-1 flows into the medial branch G7-2. With increasing Re, the velocity profile in G6-1 skews significantly towards to the medial branch G7-2. Since the flow is weak in lateral branch G7-1, the contracted tube makes the velocity profile almost symmetric after the ‘‘throat’’ which can be seen from the enlarged view in Fig. 4(c). The recirculation exists at both the entrance of G8-1 and G8-2, even the re-circulation is very weak. The path lines in
Fig. 4. The flow patterns of model 3: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) enlarged view; (d) path lines of Re ¼ 1500.
ARTICLE IN PRESS 2748
X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
Fig. 4(d) clearly show that the obstruction in G7-1 not only generates a pair of almost symmetric re-circulation in daughter airways, but also results in a stagnation region in upstream G6-1. This indicates that, in bifurcating airways, the obstruction affects the flow in both down- and up-stream airways. Fig. 5 shows the flow patterns in the bifurcation plane for model 4 in which the G7-2 airway is contracted. With increasing Re, the flow vectors are skewed significantly towards the inner wall. The obstructed airway functions like a ‘‘jet’’, for low momentum jet flow (Re ¼ 300), the velocity profile is almost symmetric at the outlets G8-3 and G8-4; for high momentum jet flow (Re ¼ 1500), the velocity profile is skewed towards the inner wall, indicating that the flow does not have enough length to recover due to the high momentum. Because of the resistance of the obstructed airway, most of the fluid
in G6-1 flows into the lateral branches G7-1, even the velocity profile skews towards the medial branch G7-2 in the second generation branch G6-1. From the path line in Fig. 5(d), most of the jet fluid flows into G8-3, and the obstruction generates the re-circulation both down- and up-stream airways. Particularly, the stagnation in up-stream G6-1 locates at almost the same position as that of model 3. 3.2. Mass flow partition In a health lung airway, the air should spread uniformly to each alveolus. However, in an obstructed airway, the air mass flow rate is biased at the first conjunction due to the non-uniform resistance downstream. The imbalance of the mass flow affects the respiration extremely. Therefore, the mass flow partition
Fig. 5. The flow patterns of model 4: (a) flow vectors of Re ¼ 300; (b) flow vectors of Re ¼ 1500; (c) enlarged view; (d) path lines of Re ¼ 1500.
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
is crucial to understand the effect of COPD on the respiratory flow. In this study, the ratio of mass flow rate of branches i to j is defined as _i m , _j m
0.8
(1)
where m ¯ i;j is the ratio of mass flow rate of branches i to j, i and j indicate the branch number, e.g. 5, 6-1, 6-2, 7-1 _ i and m _ j are the mass flow rates through branches etc., m i and j, respectively. Fig. 6 shows the influence of the obstruction on the flow partition in the obstructed branches, i.e. the flow rate in G6-1 to the total flow rate in G5. In symmetric airway (model 1), the air flow rate in G6-1 takes half of the total flow rate, indicating a balanced or symmetric flow rate between the two big branches, G6-1 and G6-2. When the second generation G6-1 is obstructed, the flow rate in G6-1 reduces significantly, and the ratio is only about 25%. The obstruction in the third generation (G71 and G7-2) results in a reduction of flow rate in branch G6-1, and the influence of the lateral branch (G7-1) and medial branch (G7-2) is almost identical. For the air flow ratio between the third generation (G7-1) and second generation (G6-1), as shown in Fig. 7, the Re may have significant influence on the flow partition. In the symmetric airway (model 1), the lateral branch G7-1 takes slightly less air compared to the medial branch G7-2 due to the centrifugal effect. When G6-1 is obstructed (model 2), the air flow ratio decreases significantly with increasing Re. It is not surprising, because the high momentum fluid induces stronger recirculation at entrance of G7-1, which blocks the fluid from entering it. It is understandable that, a small amount of air flows into G7-1 when the third generation
0.5
m6-1,5
0.4 Model 1 Model 2 Model 3 0.3
Model 1 Model 2 Model 3 Model 4
0.6
m7-1,6 -1
m ¯ i;j ¼
2749
0.4
0.2
500
1000
1500
Re Fig. 7. The variation of mass flow ratio of m ¯ 71;61 with Re.
G7-1 is obstructed; and most of the air enters G7-1 when G7-2 is obstructed. Air flow ratio is the effective way to evaluate the influence of fluid momentum on the re-circulation blockage. Figs. 8(a) and (b) show the variation of m ¯ 81;71 and m ¯ 84;72 with Re, respectively. Both the G8-1 and G8-4 airways are lateral branches for G6-1 branch, if the velocity profile in G6-1 is symmetric, i.e., it is either uniform or parabolic, the m ¯ 81;71 and m ¯ 84;72 should be the same in the symmetric airways. However, due to the centrifugal effect in G6-1 and obstruction in airways, the m ¯ 81;71 and m ¯ 84;72 are totally different and strongly dependent on Re. When G6-1 is obstructed (model 2), with increasing Re, m ¯ 81;71 increases greatly, indicating that the re-circulation at the entrance of G8-2 becomes stronger (Refer to Fig. 3(d)) which blocks the air from entering G8-2; with increasing momentum, the m ¯ 84;72 decreases significantly, the reason is the re-circulation in the mid-airway of G7-2 becomes stronger and blocks the air from getting into the G8-4 airway (Refer to Fig. 3(d)). When G7-1 is obstructed (model 3), both m ¯ 81;71 and m ¯ 84;72 increase with increasing Re. When G7-2 is obstructed (model 4), the momentum has little influence on m ¯ 81;71 , but m ¯ 84;72 drops dramatically with increasing Re.
Model 4
3.3. Pressure drop behavior
0.2
500
1000
1500
Re
Fig. 6. The variation of mass flow ratio of m ¯ 61;5 with Re.
The pressure drop in lung airways is the driving force of the respiratory flow, in the present calculation, the inlet mass flow rate is fixed and the outlet pressure is constant, therefore the obstructed airway may result in different flow resistance. The pressure drop coefficient
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
2750
Model 1 Model 2 Model 3 Model 4
Model 1 Model 2 Model 3 Model 4
3
Cp
m8 -1,7 -1
0.6
2 0.5
1 500
1000
(a)
1500
500
Re
1000
1500
Re Fig. 9. The variation of pressure drop coefficient C P with Re.
the biggest and this results in a much higher flow resistance, particularly for the low momentum fluid. This also explains why the COPD patients feel difficult to take breath.
m8 -4,7 -2
0.5
4. Conclusions 0.4
The effect of COPD on respiratory flow is numerically studied on four different four-generation models. The fully three-dimensional incompressible laminar Navier– Stokes equations are solved by a CFD solver. The numerical calculation leads to the following conclusions:
Model 1 Model 2 Model 3 Model 4 0.3 500
(b)
1000
1500
Re
Fig. 8. The variation of mass flow ratio of m ¯ 81;71 and m ¯ 84;72 with Re: (a) m ¯ 81;71 ; (b) m ¯ 84;72 .
C P is defined as CP ¼
Dptotal pin;dynamic
,
(2)
where Dptotal is the mass-weighted average total pressure drop from inlet to outlets, pin;dynamic is mass-weighted average dynamic pressure at the inlet. The variation of C P with Re for different model is plotted in Fig. 9. Generally, the pressure drop coefficient C P decreases greatly with increasing Re, it is not surprising since it is consistent with the Bernoulli’s equation. The C P distribution for models 1, 3 and 4 are quite similar, but the C P of model 2 is much higher than that of other models. The most likely explanation of this phenomenon is that the diameter reduction of model 2 is
(1) The obstructed airway alters the flow field significantly, a strong separation region exists behind the ‘‘throat’’ at higher Reynolds number. In a bifurcation airway, the obstruction may generate recirculation both upstream and downstream. These re-circulation cells block the air from entering the downstream branches which indicates that a deep breath may block the oxygen from reaching the alveoli for a COPD patient. (2) The obstructed airway has significant influence on the daughter branches and it increases the flow resistance significantly.
Acknowledgements Support given by the Research Grants Council of the Government of the HKSAR under Grant No. PolyU 5273/04E PolyU and by the Hong Kong Polytechnic University under Central Research Grant Nos. A-PG08 and G-T677 is gratefully acknowledged.
ARTICLE IN PRESS X.L. Yang et al. / Journal of Biomechanics 39 (2006) 2743–2751
References Balashazy, I., Hofmann, W., Heistracher, T., 1999. Computation of local enhancement factors for the quantification of particle deposition patterns in airway bifurcations. Journal of Aerosol Science 30, 185–203. Comer, J.K., Kleinstreuer, C., Zhang, Z., 2001. Flow structures and particle deposition patterns in double-bifurcation airway models, Part 1, air flow fileds. Journal of Fluid Mechanics 435, 25–54. Darquenne, C., 2001. Aerosol transport in the lung. In: Prisk, G.K., Paiva, M., West, J.B. (Eds.), Gravity and the Lung: Lessons from Microgravity. Marcel Dekker, Inc., New York, pp. 117–148. Hegedus, C.J., Balashazy, I., Farkas, A., 2004. Detailed mathematical description of the geometry of airway bifurcations. Respiratory Physiology & Neurobiology 141, 99–114. Isabey, D., Chang, H.K., 1981. Steady and unsteady pressure–flow relationships in central airways. Journal of Applied Physiology 51, 1338–1348. Lee, J.W., Lee, D.Y., Kim, W.S., 2000. Dispersion an aerosol bolus in a double bifurcation. Journal of Aerosol Science 31, 491–505. Liu, Y., So, R.M.C., Zhang, C.H., 2002. Modeling the bifurcating flow in a human lung airway. Journal of Biomechanics 35, 465–473.
2751
Liu, Y., So, R.M.C., Zhang, C.H., 2003. Modeling the bifurcating flow in an asymmetric human lung airway. Journal of Biomechanics 35, 951–959. Martonen, T.B., Guan, X., Schreck, R.M., 1994. Effects of carinal ridge shapes on lung airstreams. Aerosol Science and Technology 21, 119–136. Miguel, A.F., Reis, A.H., Aydin, M., 2004. Aerosol particle deposition and distribution in bifurcating ventilation ducts. Journal of Hazardous Materials B 116, 249–255. Pedley, T.J., Schroter, R.C., Sudlow, M.F., 1977. Gas flow and mixing in the airways. In: West, J. (Ed.), Bioengineering Aspects of the Lung. Marcel Dekker, NewYork, pp. 163–265. Schroter, R.C., Sudlow, M.F., 1969. Flow patterns in models of the human bronchial airways. Respiratory Physiology 7, 341–355. Snyder, B., Olson, D.E., 1989. Flow development in a model airway bronchus. Journal of Fluid Mechanics 207, 378–392. Weibel, E.R., 1963. Morphometry of the Human Lung. Springer, Academic Press, Berlin, New York. Zhao, Y., Lieber, B.B., 1994. Steady inspiratory flow in a model symmetric bifurcation. ASME Journal of Biomechanical Engineering 116, 488–496.