Table of Contents Author Guidelines Submit a Manuscript
Modelling and Simulation in Engineering
Volume 2011, Article ID 510472, 9 pages
Research Article

Simulation of Pharyngeal Airway Interaction with Air Flow Using Low-Re Turbulence Model

School of Aerospace, Mechanical and Manufacturing Engineering, RMIT University, Bundoora, VIC 3083, Australia

Received 15 October 2010; Accepted 14 February 2011

Academic Editor: Guan Yeoh

Copyright © 2011 M. R. Rasani et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


This paper aims to simulate the interaction between a simplified tongue replica with expiratory air flow considering the flow in the pharyngeal airway to be turbulent. A three-dimensional model with a low-Re SST turbulence model is adopted. An Arbitrary Eulerian-Lagrangian description for the fluid governing equation is coupled with the Lagrangian structural solver via a partitioned approach, allowing deformation of the fluid domain to be captured. Both the three-dimensional flow features and collapsibility of the tongue are presented. In addition, examining initial constriction height ranging from 0.8 mm to 11.0 mm and tongue replica modulus from 1.25 MPa to 2.25 MPa, the influence of both of these parameters on the flow rate and collapsibility of the tongue is also investigated and discussed. Numerical simulations confirm expected predisposition of apneic patients with narrower airway opening to flow obstruction and suggest much severe tongue collapsibility if the pharyngeal flow regime is turbulent compared to laminar.

1. Introduction

Collapse of the human pharyngeal airways during sleep has serious health complications with an estimated 10% of snorers being at risk to obstructive sleep apnoea [1, 2]. These episodes of partial or full cessation of breathing per hour, influence the quality of sleep, reduce brain oxygen saturation and have been linked to hypertension and heart failures [1]. The physiological mechanisms of such conditions are very closely related to flow in compliant tubes or channels. Extensive experiments of flow in Starling resistors have revealed rich variety of dynamics involved in such compressed collapsible tube systems [3, 4].

Not surprisingly, many numerical models have been developed to simulate the system and theorise the mechanisms involved in the collapse and self-excited oscillations of these compliant walls seen in the experiments. Earlier models include lumped parameter models [5] which integrated the flow variables along the whole vessel and various one-dimensional models (e.g., [68]) taking into account variation of flow variables in the longitudinal direction along the vessel but assume those variables do not vary within each cross-section. More complicated two-dimensional models have also been developed, such as those by Luo and Pedley [6, 7, 9, 10], revealing certain conditions (function of Reynolds Number and membrane tension), where steady solution are stable and below which self-excited oscillation and vorticity waves are generated in the flow. More recent three-dimensional modelling of collapsible tubes captured strong buckling of the tube [11] and also revealed both possible mechanism for flow-wall energy extraction and critical Reynolds Number for onset of self-excited oscillations [12].

Applying these developments of flow in collapsible tubes to actual physiological phenomena presents a real challenge [1]. On the subject of obstructive sleep apnoea, Chouly et al. [2, 13] and Van Hirtum et al. [14] developed an asymptotic Navier-Stokes equation coupled with linear elastic shell elements to model flow-induced deformation of a simplified tongue subjected to expiratory flow. The collapse of the tongue onto the pharyngeal wall was also simulated via in-vitro experiments to validate their model assumptions and numerical results. In general, narrowing of the airway promotes increased transmural pressure via a venturi effect, resulting in partial collapse of the airway and a nonlinear flow rate retardation as the intraluminal pressure difference is increased—a typical observation in collapsed channel called flow rate limitation.

The coupling between a two-dimensional, quasisteady, laminar, incompressible fluid flow with a shell model developed by Chouly et al. [2] was successful in obtaining efficient real-time results and was validated using pressurised latex in their experiments. Indeed, based on Reynolds Number, the flow is expected to be within the laminar (Re<2000) or transitional (2000<Re<10000) regime. However, the complex morphology of the pharyngeal airway, could give rise to three-dimensional flow features which triggers transition or turbulence regimes at lower Re [15]. In fact, observation by Hahn [16] suggests turbulence kinetic energy at approximately 5% of the inlet velocity is not unexpected in the nasal cavity.

Both neurological and physiological factors have been implicated with apneic syndrome. Previous studies (e.g., [17]) have shown that apneic patients have in general, a narrower opening in the oropharynx region (perhaps, due to obesity or tissue buildup). Pharyngeal compliance among apneic and nonapneic subjects have also been shown to vary significantly, especially during expiratory flow [17].

Hence, the present study intends to instead consider a three-dimensional model of the expiratory fluid flow using a low-Re turbulent model coupled to similar replication of the tongue. The effect of initial airway opening and the tongue stiffness on the collapsibility and flow inside the pharyngeal airway is investigated parametrically. Three-dimensional flow features and turbulence influence are also discussed.

2. Mathematical Model

2.1. Governing Equation for Fluid

Following the approach of Chouly and coworkers [2], a simplified three-dimensional flow configuration representing the pharyngeal airway is illustrated in Figure 1.

Figure 1: Pharyngeal airway model.

The converging flow through the narrowest opening in the oropharynx (at the base of the tongue) is expected to produce a jet stream and flow separation downstream of the constriction which is characterised by shear instability-induced turbulence [18]. Avoiding the intensive computational efforts involved with a three-dimensional Large Eddy Simulation (LES) and Direct Numerical Simulation (DNS), a more practical approach is employed. Steady, incompressible, turbulent flow in the three-dimensional tube is described using the Reynolds-Averaged Navier-Stokes (RANS) equations coupled with a Shear Stress Transport (SST) 𝑘-𝜔 turbulent model:𝜕𝑢𝑖𝜕𝑥𝑖𝜕𝑢=0,(1)𝑗𝑢𝑖𝜕𝑥𝑗1=𝜌𝜕𝑝𝜕𝑥𝑖+𝜕𝜕𝑥𝑗𝜈𝜕𝑢𝑖𝜕𝑥𝑗+𝜕𝑢𝑗𝜕𝑥𝑖+𝜕𝜏𝑅𝑖𝑗𝜕𝑥𝑗,(2) where 𝑢𝑖 and 𝑝 are the mean- or time-averaged velocities and pressure (𝑖=1,2,3 for three-dimensional analysis), 𝜌 is the fluid density and 𝜏𝑅𝑖𝑗 is the turbulent Reynolds stress tensor 𝜏𝑅𝑖𝑗=𝑢𝑖𝑢𝑗,(3) where 𝑢𝑖 represents the turbulent fluctuating velocity. Using Boussinesq’s assumption, this Reynolds stress tensor is related to the mean velocities by the turbulent eddy viscosity: 𝜏𝑅𝑖𝑗=𝑢𝑖𝑢𝑗=𝜈𝑇𝜕𝑢𝑖𝜕𝑥𝑗+𝜕𝑢𝑗𝜕𝑥𝑖23𝜕𝑢𝑘𝜕𝑥𝑘𝛿𝑖𝑗,(4) where 𝜈𝑇=𝑘/𝜔 is the turbulent eddy viscosity.

Therefore, in order to close the RANS equations (1) and (2), only the turbulent viscosity is required, which is determined by using 2 unknowns, that is, turbulent kinetic energy 𝑘 and turbulent rate of dissipation 𝜔. The principle behind SST 𝑘-𝜔 model is to blend the original 𝑘-𝜔 and 𝑘-𝜀 together such that, near the wall, the accuracy and robustness of the 𝑘-𝜔 model is captured, while away from the wall (i.e., free stream or free shear regions), the accuracy of the 𝑘-𝜀 is captured. Thus, the additional two equations required to solve for scalar variables 𝑘 and 𝜔 are [19]𝜕𝑢𝑗𝑘𝜕𝑥𝑗=𝜕𝜕𝑥𝑗𝜈𝜈+𝑇𝜎𝑘𝜕𝑘𝜕𝑥𝑗+𝑃𝑘𝜕𝑢0.09𝑘𝜔,(5)𝑗𝜔𝜕𝑥𝑗=𝜕𝜕𝑥𝑗𝜈𝜈+𝑇𝜎𝜔𝜕𝜔𝜕𝑥𝑗+21𝐹11𝜎𝜔2𝜔𝜕𝑘𝜕𝑥𝑖𝜕𝜔𝜕𝑥𝑖𝜔+𝛼𝑘𝑃𝑘𝛽𝜔2,(6) where 𝜎𝑘, 𝜎𝜔, 𝜎𝜔2, 𝛼, and 𝛽 are model constants (blended from the original model 𝑘-𝜔 and 𝑘-𝜀 constants), 𝐹1 represents the blending function which is dependant on the wall distance, 𝑃𝑘 is the production of turbulence due to viscous forces (function of the mean velocities), and, thus, the last two terms in both equations are effectively representing the production and destruction of each 𝑘 and 𝜔.

In addition, the SST model is further refined by employing a blending of the turbulent eddy viscosity formulation, allowing more accurate prediction of the eddy viscosity near and away from the walls: 𝜈𝑇=𝑎1𝑘𝑎max1𝜔,𝑆𝐹2,(7) where 𝑎1 is a model constant, 𝐹2 is another blending function (similarly depending on the wall distance) and 𝑆 represents an invariant measure of the strain rate [19].

Equations (1), (2), (5), (6), and (7) represent the complete fluid model which is solved by the fluid solver. A low-Re simulation with appropriate near-wall behaviour are readily captured in the SST 𝑘-𝜔  model.

Assuming the flow is symmetrical about its mid-vertical plane, a half-model of the problem is constructed as depicted in Figure 2. Following similar notations by Luo and Pedley [6], the conditions imposed at each boundary are summarised as follows:Inlet: 𝑝=𝑝in, 𝑥=(𝐿𝑢+𝐿/2),Outlet: 𝑝=𝑝out, 𝑥=(𝐿/2+𝐿𝑑),Rigid Side Walls: 𝑢𝑖=0, 𝑦2+𝑧2=𝑑,Rigid Bottom Wall: 𝑢𝑖=0, 𝑦=,Elastic Wall: 𝑢𝑖=𝑈𝑖, on Γ, Symmetrical Wall: 𝑢3=0, 𝑧=0.

Figure 2: Model definition (a) isometric view, (b) side view looking at symmetrical plane, (c) frontal view looking from the inlet.

No slip boundary conditions at the walls implies that flow velocities 𝑢𝑖 are zero at the rigid walls and at the elastic segment, 𝑢𝑖 must match the elastic wall velocities 𝑈𝑖. At the inlet, a uniform upstream pressure 𝑝in is specified, which implies a parabolic velocity profile definition at the inlet. While at the outlet, downstream pressure 𝑝out has been fixed to zero implying a typical stress-free condition.

2.2. Governing Equation for Structure

For a continuum undergoing steady deformation, Cauchy’s equation [20] is effectively a static equilibrium equation: 𝜎𝑖𝑗+𝐟=0,(8) where 𝜎𝑖𝑗 is the stress tensor and 𝐟 is the external forcing term (𝑖,𝑗=1,2,3 for three-dimensional structures) which in this case, includes external pressure 𝑝𝑒 plus the fluid pressure and fluid shear at the interface Γ. Mechanical properties of the tongue is not homogenous and anisotropic, more so with varying degrees of muscle activation. For simplicity, a homogenous isotropic material is assumed for the tongue, thus, 𝝈=𝐃𝜺 where 𝜺 is the strain vector and 𝐃 is the constitutive elastic stiffness matrix (which is a function of Young’s modulus, 𝐸, and Poisson ratio, 𝜈 of the structural material).

Considering small deformation theory, second-order terms in the Green-Lagrange strain tensor [20] could be neglected. Thus, the strain-displacement relationship could be described as: 𝜀𝑖𝑗=12𝜕𝑑𝑖𝜕𝑥𝑗+𝜕𝑑𝑗𝜕𝑥𝑖,(9) where 𝑑𝑖 is the displacement tensor. Effectively, (8) and (9) represent the solid model, capturing the elastic wall deformation under external loads.

The tongue anatomy mainly consists of water (with reported density of 1040 kg/m3 [21]) and hence, is generally accepted as incompressible [22]. Accounting for that incompressibility, a Poisson ratio 𝜈 of 0.499 has been used [2, 22] and is also adopted here. Tissues along pharyngeal airway have been reported with a range of moduli: 10–30 kPa for the vocal folds [23], 12–25 kPa for the soft palate [24] and 6 kPa [25] to 15 kPa [22] have been estimated for the tongue. In order to replicate the response of a real tongue, a hydrostatically pressurized latex shell tube was utilized by Chouly et al. [2, 13]. Similarly, a simulated shell modulus 𝐸 in combination with an external pressure 𝑝𝑒 that mimics this response is followed in this paper.

The elastic wall is assumed fixed everywhere except at the common face interfacing with the elastic segment of the fluid domain where an external pressure 𝑝𝑒 is imposed simultaneously with the pressure 𝑝 applied from the fluid domain.

3. Computational Method

A commercial finite volume solver (CFX) is used to solve the RANS equations in (1), (2), (5), (6), and (7). The fluid domain is subdivided into 8-noded three-dimensional hexahedral elements leading to an assembly of discretised algebraic equations in terms of the nodal unknowns 𝑢𝑖, 𝑝, 𝑘, and 𝜔. In general, finer meshes are located closer to the walls (elastic, rigid side and bottom) where higher velocity gradients are expected and flow separation or reattachment normally occur. In addition, denser meshes are also located in the vicinity of the constriction where sharp changes in pressure are expected.

The governing fluid equations are based on a fixed Eulerian frame of reference in space. In order to account for boundary deformation and thus, deformation of the fluid mesh, an Arbitrary Lagrangian-Eulerian (ALE) description is employed in CFX. Effectively, in an ALE framework, the mesh velocity ̃𝑢𝑗 is subtracted from the material velocity 𝑢𝑗 in the differential operator [20, 26] on the left-hand side of transport (2), (5), and (6). Thus, an additional nodal unknown (i.e., mesh velocity ̃𝑢𝑗) needs to be solved. CFX employs a diffusion model, which smoothes the internal domain nodal displacement according to a Laplace differential equation [19]: Γdisp𝛿=0,(10) where 𝛿 is the nodal displacement relative to previous mesh positions, Γdisp is some mesh stiffening parameter and (10) is solved subjected to specified nodal movement at the boundaries.

The commercial finite element solver (ANSYS) is used to solve the partial differential equations (8) and (9). Similarly, the structural domain is discretized into 8-noded SOLID185 elements and (9) is approximated using interpolation functions in terms of the unknown nodal deformations. Minimization of the variational total potential energy or weighted Galerkin residuals lead to an assembly of algebraic equation which could be solved as a function of imposed boundary conditions.

The fluid-structure interaction is achieved by satisfying either a velocity or displacement continuity (i.e., (11) or (12), resp.) and force equilibrium (13), at the common interface between both fluid and solid domain:𝑢𝑖=𝑈𝑖𝑠onΓ,(11)𝑖=𝑑𝑖onΓ,(12)𝑝+𝜏𝑖𝑗=𝜎𝑖𝑗onΓ,(13) where 𝑢𝑖 and 𝑠𝑖 are, respectively, the flow velocity and boundary displacement on the interface, 𝑈𝑖 and 𝑑𝑖 are the structure velocity and displacement, respectively, on the interface. The terms on the left-hand side of (13) represent the fluid state of stress and the right-hand side is the stress on the structure.

In order to match these conditions simultaneously in both fluid and solid solvers, a successive iteration is employed in an ANSYS Workbench platform. The fluid variables are solved based on the initial geometrical configuration. The fluid pressure and shear forces at the elastic segment are interpolated to the nodes on the elastic wall. Applying them together with the external pressure and imposed boundary conditions, the solid solver solves for the deformation of the elastic wall. These nodal deflections are then interpolated back to nodes on the elastic boundary of the fluid, which is used to effect the mesh deformation on the fluid domain. The fluid solver then solves the unknown fluid variables using the current geometrical configuration and the process is repeated until the nodal deflection and forces in (12) and (13) from current and previous coupling iterations are within a specified tolerance. The overall process can be summarised as in Figure 3.

Figure 3: Fluid-structure coupling flowchart.

Note that for steady conditions (where inertial effects are neglected and structural deformation reaches static equilibrium), elastic wall velocities 𝑈𝑖 approaches zero.

Several runs with different downstream length 𝐿𝑑 suggest that numerical results are not sensitive to the outlet proximity. In addition, finer mesh discretization also indicates similar overall wall shear distribution and more importantly, similar separation location off the elastic wall, suggesting adequate discretization.

4. Results and Discussion

In order to validate the structural model, a Poisson ratio of 0.499 and a Young’s modulus 𝐸=1.75 MPa for the isotropic elastic wall was used, giving similar load-deflection response of the tongue replica to external pressure 𝑝𝑒, as expected in [2] (Figure 4).

Figure 4: Pressure versus deflection curve for elastic wall. (Deflection is measured as lowest point on the symmetric plane. Normalization constants 𝑃𝑜=1000 Pa and 𝑜=1.6 mm.)
4.1. Parametric Investigation

For the purpose of studying the influence of the following physiological variations on the collapsibility and flow pattern in the pharyngeal airway, the following range of parameters were simulated:(i)initial constriction height ()0.8-11.0 mm,(ii)elastic wall modulus (𝐸)1.25-2.25 MPa.

The initial constriction height is varied to simulate the effect of different degrees of opening behind the tongue within the oropharynx in apneic or nonapneic subjects. The opening gap behind the tongue is typically 11 mm [27] but is perhaps lower when subjects are in sleeping position. Preliminary study by Chouly et al. [18] reported an opening of between 1 to 2 mm during sleep supine position. In line with this, variations between 0.8 to 11.0 mm have been investigated in this paper.

The influence of the modulus of the tongue is also of interest. Variation of mechanical properties in the human population is not unexpected. The degree of variation in pharyngeal compliance measured by Brown et al. [17] among apneic patients could offer some insight into the variation of stiffness in the retroglossal region. Variations in the tongue stiffness is also indicative of the degree of muscle activation [21] and have been used to capture the directional activation of muscle tissues. In this paper, considering the replication of the tongue by an elastic shell, a±0.5 MPa variation is considered.

In order to investigate the first parameter, the elastic modulus 𝐸 and external pressure 𝑝𝑒 was fixed to 1.75 MPa and 200 Pa, respectively, while the initial constriction height 𝑜 was varied. The nonlinear variation of the flow rate with intraluminal pressure difference Δ𝑝=𝑝in𝑝out is shown in Figure 5(a) for various 𝑜, revealing typical flow rate limitation phenomena with increased expiratory efforts. In particular, the increased flow rate and reduced degree of flow plateauing, suggest much reduced flow obstruction or collapse with increasing 𝑜. In addition, the lower degree of flow plateauing for the 𝑜=11 mm case would suggest a near-rigid flow rate is achieved and demonstrates the relatively low obstruction in nonapneic patients.

Figure 5: Parametric investigation (a) effect of initial airway opening 𝑜 on flow rate 𝑄 and constriction height for a fixed tongue replica modulus 𝐸=1.75 MPa, (b) effect of tongue modulus 𝐸 on flow rate 𝑄 and constriction height for a fixed initial opening 𝑜=1.2 mm (𝑄𝑜=10 L/min, 𝑃𝑜=100 Pa, and 𝑜=initialconstrictionheight).

The second parameter investigated is the effect of tongue modulus on the expiratory flow interaction with the tongue. In order to investigate this parameter, a configuration with initial constriction height of 1.2 mm and 𝑝𝑒=200 Pa was used. The effect of softening and stiffening of the tongue modulus is summarized in Figure 5(b). Reduction in modulus translates to increased compliance of the airway cross-section to the applied pressures (as suggested by reduced constriction height in Figure 5(b)). Thus, exarcebating its susceptibility to collapse and increases resistance to flow rate. The approximately linear trend of collapse with increasing Δ𝑝 is consistent with the venturi effect, where the reduction in flow area underneath the elastic wall increases the flow velocity, thus, generating much lower mean pressure (suction) that further collapses the wall.

Comparing the relative influence of both parameters, 𝑜 and 𝐸, to the wall collapse and flow rate, suggest that both parameters are important. For instance, a 40% stiffening of the tongue modulus (𝐸=1.25 to 1.75 MPa) leads to approximately 11% reduction in wall collapse, as opposed to 10% reduction in collapse for a 50% increase in initial airway opening (𝑜=0.8 to 1.2 mm).

Based on the hydraulic diameter, 𝐷=4𝐴/𝑃 where 𝐴 is taken as the cross-sectional area at the inlet and 𝑃 is the inlet perimeter, the Reynolds number, Re simulated in this investigation ranges from 274–9520. Considering a transitional flow regime in the pharyngeal airway replica (by using a low-Re turbulent model), similar to those observed in previous laminar investigations, narrower opening along the airway predisposes patients to larger wall collapse and hence, more severe flow rate limitation.

4.2. Elastic Wall Deflection

The cross-section along the narrowest opening (for a 𝑜=1.2 mm case) is shown in Figure 6 for several intraluminal pressure differences Δ𝑃. As expected, in addition to increased deflection as Δ𝑃 is increased, the profile also reveals a larger area close to the side walls where the tongue replica is supported and tends to narrow nonlinearly towards the symmetric plane 𝑧=0. This profile perhaps influences the flow patterns and will be discussed in the next section. The difference in 𝑦-position at the side wall (𝑧=0.0125) is due to the different downstream location of narrowest opening. This is associated to the downstream deflection of the elastic wall as pressure is applied in the axial direction to the upstream portion of the wall.

Figure 6: Deflection profile in the lateral direction across narrowest opening (see Figure 2(c)).

In comparison to laminar flow, the low-Re turbulent flow predicts a more severe elastic wall collapse, as illustrated in Figure 7. This would suggest a much severe pressure drop is experienced by the elastic wall in transitional or turbulent flow as oppose to laminar flow. Indeed, examining the pressure contour in Figure 8, minimum suction pressure of −81 Pa underneath the constriction is predicted in the low-Re turbulent flow compared to approximately −30 Pa in the laminar flow. This seems to be consistent with numerical findings on a rigid pharynx by Shome et al. [15] where at similar flow rates, a 40% increase in pressure drop is predicted in turbulent flow compared to if the flow is laminar and this is expected to be more severe if the flow boundaries are collapsing, as indicated here. Thus, flow regime in the pharyngeal airway has significant effect on the degree of apneic obstruction and is important to be considered in pre-operative treatments.

Figure 7: Elastic wall deflection with increasing Δ𝑃 for a 𝑜=1.2 mm case. (Normalization constants: 𝑃𝑜=100 Pa and 𝑜=initialconstrictionheight.)
Figure 8: Pressure distribution in the vicinity of the constriction at Δ𝑃=200 Pa for case with initial constriction height = 1.2 mm and 𝑝𝑒=200 Pa (𝑃𝑜=100 Pa, 𝜈𝑜=10 m/s and 𝑜=15 mm [2]).
4.3. Three-Dimension Flow Features

Typical flow patterns in the idealised airway are captured in Figure 9. In general, flow separates from the elastic wall, the bottom wall and the side walls. Three-dimensional recirculation is observed immediately downstream of the elastic wall above the jet streams in all Δ𝑃. Not surprisingly, swirling about the 𝑥-axis is also formed as the jet separates off the rigid bottom wall and is forced to follow the circular cross-section of the idealised airway.

Figure 9: Streamline plot and contours of velocity along vertical symmetry plane, bottom wall and several axial cross-sections for (a) Δ𝑃=20 Pa and (b) Δ𝑃=200 Pa (𝑜=1.2 mm case).

Streamlines from upper and bottom half of the inlet reveals the tendency for the stream to migrate towards the side walls, where the cross-sectional area is more open. This side stream (as is also visibly shown by the velocity contours peaking towards the side) is more pronounced at higher Δ𝑃 where the elastic wall deflection is further skewed downstream. Thus, the swirling strength of the inner recirculation core is increased as the outer streams follow the curvature of the walls.

Consequently, the pressure profile is also influenced by the lateral distribution of the cross-sectional opening and velocity. This is evident from the pressure profile plots along axial planes located at several lateral locations shown in Figure 10 where the critical suction pressure in the vicinity of the constriction is somewhat reduced closer to the side walls where a larger opening is anticipated.

Figure 10: Pressure profile along several axial planes at top wall.

The influence of the side wall onto the lateral pressure and velocity profile is also expected to influence the location of flow separation off both the tongue replica (elastic wall) and oropharynx base (rigid bottom wall). Figure 11 presents the axial shear profile along the bottom walls at several lateral locations 𝑧=0, 𝑧=0.005 and 𝑧=0.0, revealing a shift in the position of zero wall shear as the flow approaches the side walls. Thus, suggesting a delay in flow separation as the side wall is approached.

Figure 11: Wall shear profile along several axial planes at the bottom wall (note that elastic wall is located within 0.02136m<𝑥<0.02136m).

5. Conclusion

A three-dimensional low-Re turbulent rendering of the airway-tongue replica investigation by Chouly and coworkers [2, 13] and Van Hirtum et al. [14] have been conducted. The low-Re turbulent fluid-structure interaction simulation predicts a more severe pressure drop in the constriction compared to a laminar fluid-structure interaction, reflecting a similar trend as the turbulent rigid pharynx simulation by Shome et al. [15]. Similar flow rate limitation is observed with laminar flow, although more pronounced narrowing of the area behind the tongue is generated.

Parametric study representing some geometrical and mechanical factors that affects apnea in the airway was investigated. The investigation suggest that both a narrower opening and an increased compliance of the airway increases the susceptibility of apneic episodes. In addition, three-dimensional flow features relevant to the airway-tongue replica is also presented.

Brown et al. [17] had successfully measured changes in cross-sectional area along the pharyngeal airway during expiratory and inspiratory breathing pressures, on both apneic and nonapneic patient using an acoustic reflection technique, allowing pharyngeal compliance to be estimated as a function of the distance from the mouth. For future work, this varying stiffness along the pharyngeal airway could be incorporated in order to capture the interaction between multiple segments of differing compliance within the pharyngeal airway and how it influence collapse and reopening of the lumen. Ultimately, a real pharyngeal airway model could be used for the fluid-structure investigation.


The authors acknowledge financial support from the Government of Malaysia via the sponsorship by the Ministry of Higher Education under the IPTA Academic Training Scheme awarded to the first author. The financial support provided by the Australian Research Council (project ID LP0989452) and by RMIT University through an Emerging Researcher grant is also gratefully acknowledged. In addition, many thanks to Dr. Franz Chouly for many fruitful communications.


  1. C. D. Bertram, “Flow-induced oscillation of collapsed tubes and airway structures,” Respiratory Physiology and Neurobiology, vol. 163, no. 17#8211;3, pp. 256–265, 2008. View at Publisher · View at Google Scholar · View at Scopus
  2. F. Chouly, A. Van Hirtum, P. Y. Lagrée, X. Pelorson, and Y. Payan, “Numerical and experimental study of expiratory flow in the case of major upper airway obstructions with fluid-structure interaction,” Journal of Fluids and Structures, vol. 24, no. 2, pp. 250–269, 2008. View at Publisher · View at Google Scholar · View at Scopus
  3. C. D. Bertram, C. J. Raymond, and T. J. Pedley, “Mapping of instabilities for flow through collapsed tubes of differing length,” Journal of Fluids and Structures, vol. 4, no. 2, pp. 125–153, 1990. View at Google Scholar · View at Scopus
  4. C. D. Bertram, C. J. Raymond, and T. J. Pedley, “Application of nonlinear dynamics concepts to the analysis of self-excited oscillations of a collapsible tube conveying a fluid,” Journal of Fluids and Structures, vol. 5, no. 4, pp. 391–426, 1991. View at Google Scholar · View at Scopus
  5. A. I. Katz, Y. Chen, and A. H. Moreno, “Flow through a collapsible tube. Experimental analysis and mathematical model,” Biophysical Journal, vol. 9, no. 10, pp. 1261–1279, 1969. View at Google Scholar · View at Scopus
  6. X. Y. Luo and T. J. Pedley, “A numerical simulation of steady flow in a 2-D collapsible channel,” Journal of Fluids and Structures, vol. 9, no. 2, pp. 149–174, 1995. View at Publisher · View at Google Scholar · View at Scopus
  7. T. J. Pedley, “Longitudinal tension variation in collapsible channels: a new mechanism for the breakdown of steady flow,” Journal of Biomechanical Engineering, vol. 114, no. 1, pp. 60–67, 1992. View at Google Scholar · View at Scopus
  8. A. H. Shapiro, “Steady flow in collapsible tubes,” ASME Journal of Biomechanical Engineering, vol. 99, no. 3, pp. 126–147, 1977. View at Google Scholar · View at Scopus
  9. X. Y. Luo and T. J. Pedley, “A numerical simulation of unsteady flow in a two-dimensional collapsible channel,” Journal of Fluid Mechanics, vol. 314, pp. 191–225, 1996. View at Google Scholar · View at Scopus
  10. X. Y. Luo and T. J. Pedley, “The effects of wall inertia on flow in a two-dimensional collapsible channel,” Journal of Fluid Mechanics, vol. 363, pp. 253–280, 1998. View at Google Scholar · View at Scopus
  11. A. L. Hazel and M. Heil, “Steady finite-Reynolds-number flows in three-dimensional collapsible tubes,” Journal of Fluid Mechanics, no. 486, pp. 79–103, 2003. View at Publisher · View at Google Scholar · View at Scopus
  12. M. Heil and S. L. Waters, “How rapidly oscillating collapsible tubes extract energy from a viscous mean flow,” Journal of Fluid Mechanics, vol. 601, pp. 199–227, 2008. View at Publisher · View at Google Scholar · View at Scopus
  13. F. Chouly, A. Van Hirtum, P. Y. Lagrée, X. Pelorson, and Y. Payan, “Modelling the human pharyngeal airway: validation of numerical simulations using in vitro experiments,” Medical and Biological Engineering and Computing, vol. 47, no. 1, pp. 49–58, 2009. View at Publisher · View at Google Scholar · View at Scopus
  14. A. Van Hirtum, X. Pelorson, and P. Y. Lagrée, “In vitro validation of some flow assumptions for the prediction of the pressure distribution during obstructive sleep apnoea,” Medical and Biological Engineering and Computing, vol. 43, no. 1, pp. 162–171, 2005. View at Publisher · View at Google Scholar · View at Scopus
  15. B. Shome, L. P. Wang, M. H. Santare, A. K. Prasad, A. Z. Szeri, and D. Roberts, “Modeling of airflow in the pharynx with application to sleep Apnea,” Journal of Biomechanical Engineering, vol. 120, no. 3, pp. 416–422, 1998. View at Google Scholar · View at Scopus
  16. I. Hahn, Modelling Nasal Airflow and Olfactory Mass Transport, University of Pennsylvania, Philadelphia, Pa, USA, 1992.
  17. I. G. Brown, T. D. Bradley, and E. A. Phillipson, “Pharyngeal compliance in snoring subjects with and without obstructive sleep Apnea,” American Review of Respiratory Disease, vol. 132, no. 2, pp. 211–215, 1985. View at Google Scholar · View at Scopus
  18. F. Chouly, A. van Hirtum, P.-Y. Lagrée et al., “Simulation of the retroglossal fluid-structure interaction during obstructive sleep Apnea,” in Lecture Notes in Computer Science, M. Harders and G. Szekely, Eds., pp. 48–57, Springer, Berlin, Germany, 2006. View at Google Scholar
  19. ANSYS, ANSYS CFX Solver Theory Guide—Release 12.
  20. G. Xia and C. L. Lin, “An unstructured finite volume approach for structural dynamics in response to fluid motions,” Computers and Structures, vol. 86, no. 7-8, pp. 684–701, 2008. View at Publisher · View at Google Scholar · View at Scopus
  21. F. A. Duck, Physical Properties of Tissue: A Comprehensive Reference Book, Academic Press, London, UK, 1990.
  22. Y. Payan, G. Bettega, and B. Raphael, “A biomechanical model of the human tongue and its clinical implications,” in Lecture Notes in Computer Science, W. M. Wells, A. Colchester, and S. Delp, Eds., pp. 688–695, Springer, Berlin, Germany, 1998. View at Google Scholar
  23. Y. Min, I. Titze, and F. Alipour, “Stress-strain response of the human vocal ligament,” in NCVS Status and Progress Report, pp. 131–137, 1994. View at Google Scholar
  24. Y. Huang, D. P. White, and A. Malhotra, “The impact of anatomic manipulations on pharyngeal collapse: results from a computational model of the normal human upper airway,” Chest, vol. 128, no. 3, pp. 1324–1330, 2005. View at Publisher · View at Google Scholar · View at Scopus
  25. Y. Huang, A. Malhotra, and D. P. White, “Computational simulation of human upper airway collapse using a pressure-/state-dependent model of genioglossal muscle contraction under laminar flow conditions,” Journal of Applied Physiology, vol. 99, no. 3, pp. 1138–1148, 2005. View at Publisher · View at Google Scholar · View at Scopus
  26. J. Donea, A. Huerta, J.-P. Ponthot, and A. Rodríguez-Ferran, “Arbitrary lagrangian-eulerian methods in encyclopedia of computational mechanics,” in Fundamentals, E. Stein, R. de Borst, and T. J. R. Hughes, Eds., vol. 1, John Wiley & Sons, 2004. View at Google Scholar
  27. E. Sforza, W. Bacon, T. Weiss, A. Thibault, C. Petiau, and J. Krieger, “Upper airway collapsibility and cephalometric variables in patients with obstructive sleep Apnea,” American Journal of Respiratory and Critical Care Medicine, vol. 161, no. 2, part I, pp. 347–352, 2000. View at Google Scholar · View at Scopus