CFD simulations of respiratory airflow in human upper airways response to walking and running for oral breathing condition.

Walking and running are common types of physical activities that people do in day to day living, to improve health and physical fitness or for recreation. During a physical activity, rate and depth of breathing increase because working muscles need extra oxygen in order to produce energy. In this study, computational fluid dynamics (CFD) simulations were used to investigate respiratory airflow dynamics in human upper airways response to walking and running for oral breathing. The numerical simulations were done in a realistic CT-scan airway model using ANAYS Fluent 19.0 software. Flow fields were analysed numerically and flow patterns were investigated in the airway model during inspiration and expiration response to walking and running. The axial velocity distributions and secondary flow patterns for the two respiratory phases were analysed response to the two physical activities at different cross-sections of the airway model. The maximum velocity, wall pressure and wall shear stress values for running were respectively 3.2, 9.4 and 5.9 times higher than that of walking during inspiration. The mixing of flow streamlines was observed to be higher during running than walking because of more significant turbulence. More skewed flows at airway curvatures were observed at inspiration than expiration. The results of this study supported the fact that running is a more intense activity than walking from a respiratory dynamics point of view.

Introduction

Colditz et al. [1] mentioned that engaging in regular physical activities is one of the best ways to improve overall health. They added that physical activities could prolong life, lower the risk of numerous diseases, promote self-confidence and manage weight. During physical activity, body muscles use more energy. In order to produce this required energy, the working muscles need extra oxygen to release this energy from glucose. Hence lungs have to work harder to take fresh oxygen into the body and get rid of carbon dioxide, the waste product created when energy is produced, from the blood. So, physical activity makes breathing rate higher than usual. Walking and running are usual physical activities that people perform. Wilkin et al. [2] compared energy expenditure between walking and running. No adequate studies are reporting the respiratory airflow dynamics response to these physical activities.

The advancement of CFD tools has attracted many researchers’ attention to study airflow dynamics in a model of human airways. CFD simulations of respiratory airflow have been done in simplified airway models (symmetric or asymmetric) [3, 4] and realistic CT-based airway models [5]. Due to the complexity of the human airway structure, computing time and computer memory, most CFD simulations are limited to small parts of the respiratory tract [6]. Numerous studies of respiratory airflow have been carried out in human upper airways based on CT-scanned airway models. Wang et al. [7] used a combination of CT-scan model (consists of nasal cavity, pharynx and larynx) and Weibel model (containing the tracheobronchial tree up to third generation) for the study of airflow characteristics in the human respiratory tract. Yu et al. [8] studied the effects of airways configuration on respiratory flow in a CT-scan airway model including the nasal cavity, oral cavity, pharynx, larynx, trachea and some bronchi. To investigate inspiratory airflow, Shang et al. [9] employed a realistic respiratory airway model from nasal and oral openings to terminal bronchi up to fifteenth generations. Kou et al. [10] analysed the airflow dynamics during coughing in a CT-scan airway model from oral cavity to bronchi up to the fifth generations. Calmet et al. [11] used a CT-scan airway model from nasal cavity to bronchi up to third generations to investigate airflow dynamics during a rapid breath.

It is important to understand airflow characteristic (velocity, wall pressure and wall shear stress) and airflow patterns in human airways while breathing during high levels of physical activity. The flow is more time-dependent at high breathing frequencies during a rapid breath response to a higher level of physical activity, and transient simulations should be adopted to describe airflow dynamics in the airways [12]. During a higher level of physical activity, high flow rate and laryngeal jet create turbulent flow in human upper airways. The turbulent effects of airflow should also be considered for better description of airflow characteristics for respiratory airflow during higher levels of physical activity. Turbulence models are necessary because we cannot afford big enough computers to compute every scale of motion directly [13].

This study's objective was to investigate respiratory airflow dynamics response to walking and running in human upper airways using CFD simulations. Unsteady-state flow simulation and turbulence modeling were intended to be used in the numerical characterization of flow fields and visualisation of flow patterns.

Methods

Airway model

Tsega et al. [14] studied breathing patterns response to different physical activity levels only from mouth ventilation. The average height of subjects in this study was 166.0 cm. By scaling using this average height [15, 16], the average length and diameter of trachea for a subject were 85.9 mm and 15.6 mm respectively. Based on the small difference in diameter and length of the trachea (compared to the simplified model constructed by scaling) and its suitability to simulate only from mouth breathing, the CT-scan airway model used by Rahimi-Gorji et al. [17] was adopted in the study. For this CT-scan airway model, the trachea length is 92.38mm and its inlet hydraulic diameter is 16.57mm. The airway model consists of 13 outlets (7 in the left and 6 in the right bronchial tree) up to the sixth generations is shown in Figure 1a at different views. The cross-section of the mouth inlet of the model used in this study is not circular and its area is 1.1837cm2 which is different from the model used in [17].

Figure 1

(a) The CT-scan airway model at different views (b) Different cross-sections on the airway model.

Mesh generation

Unstructured tetrahedron meshes were generated on the airway model using ANSYS 19.0 Workbench Meshing tool for numerical computation. Inflation layers were generated to capture the boundary layer flows near the wall of the airways (Figure 2). A grid independence test was carried out using steady flow simulations with an inlet velocity of 2.18 m/s (average inlet velocity for inspiration in response to sitting at rest [14]) at the mouth starting from 1,423,478 cells (M1) and increasing the number cells to 2,185,129 (M2) and 3,410,756 (M3). The mean velocity of five points on the right bronchus was calculated for each simulation. The mean velocities were 1.7901, 1.8829 and 1.8946 m/s for mesh M1, M2 and M3, respectively. The mean velocities were almost the same for 2,185, 129 and 3,410,756 cells [18]. Taking the computational memory and time into account, mesh with 2,185,129 cells was selected for simulations in this study.

Figure 2

Mesh preview at selected locations.

Governing equations and boundary conditions

The air flowing in the human respiratory passageways can be regarded as homogeneous, Newtonian and incompressible fluid [19]. The Reynolds numbers of the inspiratory and expiratory airflow at the mouth inlet and larynx of the airway model response walking and running were calculated as shown in Table 1. The inspiratory and expiratory flow rates from [14] and the airway model geometry were used to obtain the Reynolds numbers. From this, it is observed that a turbulent flow occurred in the airway model. Using a Reynolds averaged Navier-Stokes (RANS)-based turbulence mode, the unsteady airflow in the respiratory tract can be described by the following time-averaged equations:

Table 1

Respiratory parameters (averages) used in the simulations (Adopted from [14]).

PhysicalActivity	Inspiration				Expiration
	Flow rate (L/s)	InletVelocity (m/s)	Re at inlet	Re at Larynx	Flow rate (L/s)	InletVelocity (m/s)	Re at inlet	Re at Larynx
Walking	0.7101	5.9994	5052	6129	0.6087	-5.1423	4330	5253
Running	2.3929	20.2151	17022	20651	2.2333	-18.8674	15887	19274

Continuity Equation:

Momentum Equation:where is the mean velocity component in the , and directions, is the mean pressure, is air density (), is dynamic viscosity of the air () for laminar flow [6] and is the turbulent viscosity. The k-ω shear stress transport (SST) turbulence model was employed in the simulation which is suitable to predict different flow regimes and more accurate and reliable [9, 10, 20]. The k-ω SST model combines the k-ω turbulence model to use near walls and standard k–ε turbulence model to use away from walls using blending functions [21]. In k-ω SST model, the transport equations of the turbulent kinetic energy (k) and the specific dissipation rate (ω) are given bywhere and is turbulence stress tensor given . The turbulent viscosity in k- ω SST model is computed as , where , is the strain rate tensor and the function is given by where d is distance to the nearest surface. The model constants , , and in Eq. (3) and Eq. (4) are computed aswhere is a constant in k- ω model and is the corresponding constant in the k- ω formulation of the k-ε mode. In Eq. (5), the blending function is given by

The constant is given by. The constants in the k- ω SST model have the following empirically derived values [22]:

From tidal volume, inspiration time [14] and cross-sectional area of mouth inlet, the average respiratory flow rates and the corresponding average velocities at the mouth inlet in response to walking and running were computed as shown in Table 1. Uniform inlet velocity was imposed at the mouth inlet using the average velocity for each breathing conditions [19, 23, 24]. The velocity at the mouth inlet was set to be negative for expiration [6]. A zero relative pressure was specified at each of the thirteen outlets and no-slip condition was applied at the wall of the airways [25, 26, 27].

Numerical methods

To solve the governing equations, Eqs. (1), (2), (3), and (4), and simulate airflow in the constructed airway model, the commercial CFD software ANSYS Fluent (version 19.0) was used. The SIMPLE algorithm was applied for the pressure-velocity coupling. The second order upwind difference scheme was applied for discretizing the convective terms. The diffusion gradients were approximated using Green-Gauss cell-based method. The time step of 0.0001s was used for the unsteady simulations which ensure that the Courant number is less than 1. The flow time employed for each simulation was 0.4s [28]. Convergence criteria for numerical results were residuals less than for continuity, momentum and turbulence equations. A 64-bit laptop PC with Intel Core i5 processor and 8GB RAM was used for computations.

Validation

To validate computational methodology employed in this study, numerical simulation of inspiratory airflow () was performed and velocity profile at a cross-section of a fourth generation bronchus of the right lung was compared with the corresponding experimental data of Zhao and Lieber [29]. The comparison is displayed in Figure 3 and shows that the simulation result is in good agreement with the experimental data. The discrepancies may be because of differences in geometries of the airway models.

Figure 3

Comparison between present computational works with experimental data.

Results and discussions

Velocity distributions and flow patterns

Figure 4 shows inspiratory and expiratory airflow velocity streamlines response to walking and running in the CT-scan airway model. As physical activity level increased from walking to running, the airflow velocity increased in the airway model for both inspiration and expiration. For both breathing conditions in the upper respiratory tract, high-speed flow was observed near the end of the oral cavity and larynx regions which agrees with the result predicted by [9]. In lower respiratory tract, the airflow velocity was generally higher at the outlet bronchi than at the main airways which accords with previous CFD simulations [6, 10, 11]. The maximum inspiratory airflow velocities reached in the airway model in response to walking and running were 12.15 and 38.96 m/s respectively. The corresponding maximum expiratory airflow velocities were 10.19 and 34.35 m/s. These maximum values were attained in the larynx region. The maximum velocity in response to running was 3.2 and 3.4 times greater than that of walking at inspiration and expiration, respectively.

Figure 4

Airflow velocity streamlines response to (a) walking (b) running.

Flow streamlines show the occurrence of re-circulatory flow at the portion of the oral cavity below the nasopharynx, oropharynx and larynx region [9]. This was the effect of narrowed paths in the regions. The re-circulation was observed prominently for walking at both respiratory phases. At the trachea, streamline patterns were relatively similar in response to walking and running for the corresponding respiratory phase. It was observed that more air entered into the right main bronchus than the left main bronchus during inspiration [30]. The airflow patterns for inspiration and expiration were different. The number of streamlines reached distal bronchi was higher for running than walking during inspiration as shown in Figure 4.

The axial velocity distribution contours and secondary flow velocity vectors on cross-sections 1-10 are shown in Figures 5 and 6 for inspiration and expiration response to walking and running. The cross-sections in the airway model are indicated in Figure 1b. The axial velocity distributions and magnitudes at the cross-sections were different for walking and running at both respiratory phases. The axial velocity distributions differences were relatively greater at the cross-sections 2, 4 and 7 during inspiration and at the cross-sections 3 and 4 during expiration. Generally, more skewed flows at airway curvatures were observed at inspiration than expiration [6, 9]. There were significant secondary flow variations on each cross-section between walking and running during inspiration and expiration. The more disturbing flow was observed at the cross-sections due to the secondary flow motion for running than walking at both respiratory phases. The study indicates that there are differences in pattern and velocity of airflow during walking and running.

Figure 5

Axial velocity distributions and secondary flow vectors on cross-sections 1-10 for inspiratory airflow response to (a) walking (b) running.

Figure 6

Axial velocity distributions and secondary flow vectors on cross-sections 1-10 for expiratory airflow response to (a) walking (b) running.

Pressure distribution

Figure 7 shows wall pressure distribution for respiratory airflow response to walking and running in the CT-scan airway model. For inspiration, relatively high wall pressure was observed at the oral cavity and decreased to the end bronchi for the two breathing conditions [10]. Wall pressure increased as the level of physical activity increased. The pressure drops in the airway model were 146.48 and 1373.93 Pa in response to walking and running respectively. The pressure drop in response to running was 9.4 times greater than that of walking. A rapid pressure drop was observed for both walking and running when air passing the oropharynx and flow resistance of air was relatively higher at the oral cavity. The flow resistance of air was somewhat higher at the oral cavity which is consistent with the study results of [8, 9, 30]. For expiration, pressure drops were 100.46 and 1034.73 Pa in response to walking and running respectively. The pressure drop in response to running was 10.3 times greater than that of walking. The airflow underwent a small pressure drop when passing from terminal bronchi to the main airways compared to that of inspiration as reported by [8].

Figure 7

Wall pressure distributions for respiration response to (a) walking (b) running.

Wall shear stress

To understand the interaction of the flowing air and the airway walls, it is important to investigate the wall shear stress. For inspiration and expiration response to walking and running, the wall shear stress distributions are displayed in Figure 8. The wall shear stress variations for walking and running correspond to airflow velocity on the airway model [6]. For both walking and running, the wall shear stress was higher at inspiration than the expiration. For inspiration, the maximum wall shear stress values reached were 5.52 and 32.45 Pa response to walking and running respectively. The corresponding maximum wall shear stress values for expiration were 3.44, and 25.20 Pa. The wall shear stress response to running was about 5.9 times greater than that of walking for inspiration and over 7.3 times greater for expiration. The study provides the quantitative increment of shear stress near walls of airways with increased level of physical activity.

Figure 8

Wall shear stress distributions for respiration response to (a) walking (b) running.

Turbulence kinetic energy

The distribution of turbulence kinetic energy in the vertical mid-plane of the airways is shown in Figure 9. The maximum turbulent kinetic energies reached in the airway model response to walking and running were 3.44 and 30.48 m2/s2 respectively. The corresponding maximum expiratory turbulent kinetic energies were 7.03 and 45.24 m2/s2. Because of greater turbulence, more mixing of streamlines was observed response to running than that of walking (Figure 4) [6, 31].

Figure 9

Turbulence kinetic energy for respiration response to (a) walking (b) running.

Conclusions

In this study, numerical simulation of respiratory airflow response to walking and running was performed in a realistic CT-scan airway model during oral breathing. Increasing physical activity level from walking to running caused a higher value of flow characteristics (airflow velocity, pressure and wall shear stress) in the entire human upper airways region. The changes were described quantitatively in the study. Flow patterns of respiratory airflow response to walking and running were different in human upper airways. Due to greater turbulence incidence, more mixing of flow streamlines was observed during running than walking. More skewed flows at airway curvatures response to walking and running were observed for inspiration than for expiration. The effect of secondary flow on the axial velocity distribution in human upper airways was higher for running than walking. The study can help to create awareness regarding respiratory airflow dynamics related to physical activity. The power of computational fluid dynamics techniques may be appreciated in simulating human physiological processes.

Declarations

Author contribution statement

Endalew Getnet Tsega: Conceived and designed the experiments; Performed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Funding statement

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Data availability statement

No data was used for the research described in the article.

Declaration of interests statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.