View source for Best Practice Advice AC7-02

{{ACHeader
|area=7
|number=02
}}
__TOC__
=Airflow in the human upper airways=
'''Application Challenge AC7-02'''&nbsp;&nbsp;&nbsp;© copyright ERCOFTAC 2020
=Best Practice Advice=
==Key Fluid Physics==
<!--
Briefly   describe  the  key   fluid   physics/flow   regimes   which   exert   an   influence   on   the  DOAPs. Ideally   this   should   draw   together   into   a   coherent   picture   the   associated   UFR   descriptions together with any important interactions which are AC specific. Mention the UFRs associated with   this   AC   that   you   have   considered   in   drafting   your   best   practice   advice.  ''Access   the Knowledge Base to find the UFRs associated with your AC''.
-->

In the present AC, experiments and simulations were conducted at a flowrate of 60 L/min through an upper airway geometry. 
At this flow conditions, the Reynolds number for air in the trachea is 4920, which is well within the turbulent regime. 
Geometric effects, such as the bent in the oropharyngeal region and the constriction at the laryngeal glottis (just upstream of the trachea, see [[#Figure25|figure 25]]
) enhance turbulence levels as the air moves from the inlet to the region of the trachea. 
Turbulent kinetic energy levels reach a peak in the shear layer formed between the high speed laryngeal jet and the surrounding (low speed) air (see [[#Figure25|figure 25]]). 
The characteristics of the laryngeal jet formation bear a resemblance to the flow through a constricted pipe, which can be classified as a free shear flow where the wall serves to confine the spreading of the jet rather than producing turbulence ([[Best Practice Advice AC7-02#Tawhai2011|Tawhai & Lin, 2011]]). 
High turbulence levels persist in the region of the first bifurcation (stations H1-H2 & J1-J2 in [[CFD Simulations AC7-02#Figure12|figure 12]](b)).

==Application Uncertainties==
<!--
List   any   uncertainties   which   make   a   high   fidelity   CFD   model   difficult   to   assemble.   Typical examples might include:
*a gas leakage between two components which is difficult to resolve on practical meshes, andeven if it is resolved, the leakage flow conditions may not be known.
*flow conditions at inlet to the AC (or indeed other boundaries) which may be complex and not precisely known.
*fine details of the geometry are imprecise.
Briefly discuss the sensitivity of the DOAP predictions to these uncertainties and their impact on the   BPA.   In   particular,   can   clear,   unequivocal   BPA   be   given   or   is   it   necessary   to   introduce appropriate caveats.
-->

The differences between measurements and simulations can result from several uncertainties involved in the tests. 
A first source of uncertainty are the inlet conditions, which are not perfectly matched between the measurements and the computations. 
In the experiments, the lung model was placed in an open liquid tank with a piston diaphragm pump attached to a linear actuator to achieve a quasi-stationary inspiratory flow. 
The stroke of the piston followed a cyclic triangular function with an adjustable falling constant slope and thus constant velocity to match different flow rates during inspiration. 
The measured mean velocity at the inlet of the model, shown in [[Test Data AC7-02#Figure7|figure 7]], is asymmetric, probably due to the action of the piston diaphragm pump. 
In the computations, instead of reproducing the measured inlet conditions, either uniform or turbulent inlet velocity profiles were prescribed. 
Due to a leakage flow between the upper and lower part of the model in the experiments, the achieved flowrate within the main bifurcation and bronchi region was about 10% lower than in the upper part of the model. 
As a result, a maximum flowrate of <math> Q_{w/g} </math> = 28.56 L/min could be achieved in the measurements. 
This value is slightly lower than <math> Q_{w/g} </math> = 31.75 L/min, which is the target value for an equivalent air flowrate of 60L/min through the model. 
Although the flow is well within the turbulent regime, the theoretical maximum Reynolds number decreases from 4921 to 4286.

==Computational Domain and Boundary Conditions==

The geometry of the extrathoracic airways must be included because turbulence is generated in this region that propagates in the first airway generations. 
Concerning the boundary conditions, the inlet velocity profile is important and thus realistic inlet conditions should be used. 
At the outlets, it is important to apply correct pressures such that the ventilation of the airway tree is physiologically realistic ([[Best Practice Advice AC7-02#Yin2010|Yin et al., 2010]]). 
In the present AC, in order to simplify the experimental setup and be able to perform the flow measurements, uniform pressures were prescribed at all outlets.

Concerning the inlet conditions for the turbulent variables in RANS calculations, the application of a turbulence intensity of 5% for the k-ω SST model at the extended inlet (10xDinlet) yielded higher turbulent kinetic energy values close to the inlet of the model compared to the mapped inlet condition (Inlet 1). 
The k-ε models were found to provide overall higher turbulence levels than the k-ω SST model, especially at the near-wall regions.

==Discretisation and Grid Resolution==

Since it is not possible to generate a structured hexahedral grid for the present geometry due to its complexity, a higher refinement ratio should be applied to avoid numerical diffusion. 
In addition to that, layers of prismatic elements should be added near the wall boundaries for a better prediction of this region, not only with regard to flow properties itself, but the flow conditions seen by the particles, i.e. mean velocity and turbulence properties. 
Airflow through the glottis constriction at the larynx, illustrated in [[#Figure25|figure 25]], bears a resemblance to flow through a constricted pipe. 
This type of flow may be classified as a free shear flow in which the wall serves to confine the spreading of the jet rather than producing turbulence. 
In this case, turbulence is most active at the interface between two free streams (high and low speed) and along the jet core and therefore, fine mesh resolution should be placed accordingly to capture strong turbulence activities. 
Use of a strict y+ = 1 condition in the generation of the near-wall mesh but extremely coarse mesh in the core region of the airway model is conceptually wrong ([[Best Practice Advice AC7-02#Tawhai2011|Tawhai & Lin, 2011]]).
Recommended values for the parameters involved in mesh generation (initial cell height, average expansion ratio, number of near-wall prism layers, average cell volume in the domain, number of computational cells etc.) can be found in Tables 4,5 (LES) and 6 (RANS).

LES were found to give similar results independent of the discretisation method used (Finite Volume or Finite Element).

<div id="Figure25"></div>
[[File:jet.png|center|thumb|500px|'''Figure 25''': Laryngeal jet formed at the glottis constriction. (a) Isosurfaces of <math> |<u>| </math> = 6.7m/s ; (b) Mean velocity contours; (c) Turbulent kinetic energy contours.]]

==Turbulence Models==


The ability of four LES subgrid-scale models was assessed by comparing their predictions to the PIV data. 
These models were the WALE, the Smagorinsky, the variational multiscale (VMS) WALE and the QR eddy-viscosity model from Verstappen. 
All the models were found to provide very similar results. 
It is concluded that the influence of the subgrid scales on the airflow in the human upper airways is small, and the choice of subgrid-scale turbulence model is not as important as in cases with higher Reynolds numbers.

In general the RANS prediction for the velocity magnitude follow trends similar to those of the measured data. 
RANS however show poor performance at locations with shear layers, recirculation and flow separation. 
Moreover, RANS simulations consistently underestimate turbulent kinetic energy levels compared to PIV and LES data. 
This is a direct consequence of discarding many of the elements of the underlying turbulence physics when solving only for the mean flow (see discussion in section [[Evaluation AC7-02|Evaluation]]). 
The k-ε models provides overall higher turbulence levels than the k-ω SST turbulence model.

Results presented in AC7-01 showed that the deposition fractions obtained with RANS models achieved a quite good agreement with LES and measurements. 
This might seem surprising given the tendency of RANS to underpredict turbulence intensities at several stations downstream of the glottis constriction. 
This apparent paradox probably relates to the fact that the dominant deposition mechanism in the upper airways is inertial impaction. 
While inertial impaction tends to be dominated by mean flow effects, turbulent dispersion still plays an important role, especially in regions where there is significant large-scale anisotropy in the turbulence. 
Hence, it is known that using only the time-averaged air velocity field for deposition studies (this is the case in RANS) leads to deposition overpredictions ([[Best Practice Advice AC7-02#Matida2004|Matida et al., 2004]]).
Still, RANS deposition predictions can be improved when they are used together with a turbulent dispersion model, as was done in AC7-01. 
When using a turbulent dispersion model, individual particles are allowed to interact successively with discrete eddies, each eddy having length, velocity and lifetime characteristic scales obtained from the primary flow calculation results. 
It is therefore important to correctly select the turbulent dispersion model and it’s parameters for accurate deposition predictions. 
LES on the other hand, since they resolve the large scale eddies, do not need a model to account for these in particle-laden flows. 
The study of [[#armenio1999|Armenio et al. (1999)]] has shown that the motion of inertial particles in low to moderate Reynolds number flows, is not influenced from the unresolved sub-grid scales in LES.

In conclusion, LES are more capable than RANS in predictions of airflow in the human upper airways, since they can account better for the physics of the turbulent flow without the need to adjust model parameters.

==Recommendations for Future Work==

The present application challenge focuses on the airflow that develops in the human upper airways. 
The airway geometry has been considered rigid. 
In reality, the lung expands and contracts and the airway walls deform during inhalation and exhalation ([[Best Practice Advice AC7-02#Mead-Hunter2013|Mead-Hunter et al., 2013]]). 
In addition, there is periodic movement of the glottal aperture during tidal breathing that regulates the respiratory airflow dynamics ([[Best Practice Advice AC7-02#Xi2018|Xi et al., 2018]]).
It is therefore important to assess the effect of wall deformation on the developed flow features inside the airways.

The higher temperature and humidity of the human body compared to the inhaled ambient air results in heat and water vapor transfer as the air is transported in the airways 
([[Best Practice Advice AC7-02#Wu2014|Wu et al., 2014]]).
During inhalation, air is heated and humidified by the airway walls until it reaches the body temperature and approximately 100% relative humidity. 
Such effects can play a role in air and inhaled aerosols transport inside the airways and should be further examined in the future studies.

==Acknowledgements==

The present application challenge is based upon work from COST Action MP1404 SimInhale "Simulation and pharmaceutical technologies for advanced patient-tailored inhaled medicines", supported by COST (European Cooperation in Science and Technology - www.cost.eu).

==References==
<div id="Adrian2011"></div>
Adrian, R.J. and Westerweel, J. 2011 
:Particle Image Velocimetry. ''Cambridge University Press, Cambridge.''
<br/>

<div id="armenio1999"></div>
Armenio, V., Piomelli, U. & Fiorotto, V. 1999
:Effect of the subgrid scales on particle motion. ''Physics of Fluids'' '''11''' (10), 3030&nbsp;&ndash;&nbsp;3042.
<br/>

<div id="Charnyi2017"></div>
Charnyi, Sergey, Heister, Timo, Olshanskii, Maxim A. & Rebholz, Leo G. 2017 
:On conservation laws of Navier-Stokes Galerkin discretizations. ''Journal of Computational Physics'' '''337''', 289&nbsp;&ndash;&nbsp;308.
<br/>

<div id="Codina2001"></div>
Codina, R. 2001 
:Pressure stability in fractional step finite element methods for incompressible flows. ''Journal of Computational Physics'' '''170''', 112&nbsp;&ndash;&nbsp;140.
<br/>

<div id="Hughes2000"></div>
Hughes, Thomas J.R., Mazzei, Luca & Jansen, Kenneth E. 2000 
:Large eddy simulation and the variational multiscale method. ''Computing and Visualization in Science'' '''3''' (1), 47&nbsp;&ndash;&nbsp;59.
<br/>

<div id="Janke2019"></div>
Janke, T., Koullapis, P., Kassinos, S.C. & Bauer, K. 2019 
:PIV measurements of the Siminhale benchmark case. ''European Journal of Pharmaceutical Sciences'' '''133''', 183&nbsp;&ndash;&nbsp;189.
<br/>

<div id="Jasak1996"></div>
Jasak, H. 1996 
:Error analysis and estimation for the finite volume method with applications to fluid flows. ''PhD thesis, Department of Mechanical Engineering, Imperial College of Science, Technology and Medicine, London, UK''.
<br/>

<div id="Jofre2014"></div>
Jofre, L., Lehmkuhl, O., Ventosa, J., Trias, F.X. & Oliva, A. 2014 
:Conservation properties of unstructured finite-volume mesh schemes for the Navier-Stokes equations. ''Numerical Heat Transfer, Part B: Fundamentals'' '''54''' (1), 53&nbsp;&ndash;&nbsp;79.
<br/>

<div id="Kim1985"></div>
Kim, J. & Moin, P. 1985 
:Application of a fractional-step method to incompressible Navier-Stokes equations. ''Journal of Computational Physics'' '''59''' (2), 308&nbsp;&ndash;&nbsp;323.
<br/>

<div id="Koullapis2018"></div>
Koullapis, P., Kassinos, S. C., Muela, J., Perez-Segarra, C., Rigola, J., Lehmkuhl, O., Cui, Y., Sommerfeld, M., Elcner, J., Jicha, M., Saveljic, I., Filipovic, N., Lizal, F. & Nicolaou, L. 2018 
:Regional aerosol deposition in the human airways: The SimInhale benchmark case and a critical assessment of ''in silico'' methods. ''European Journal of Pharmaceutical Sciences'' '''113''', 77&nbsp;&ndash;&nbsp;94.
<br/>

<div id="Koullapis2016"></div>
Koullapis, P. G., Kassinos, S.C., Bivolarova, M. P. & Melikov, A. K. 2016 
:Particle deposition in a realistic geometry of the human conducting airways: Effects of inlet velocity profile, inhalation flowrate and electrostatic charge. ''Journal of Biomechanics'' '''49''', 2201&nbsp;&ndash;&nbsp;2212.
<br/>

<div id="Lehmkuhl2019"></div>
Lehmkuhl, O., Houzeaux, G., Owen, H., Chrysokentis, G. & Rodriguez, I. 2019 
:A low-dissipation finite element scheme for scale resolving simulations of turbulent flows. ''Journal of Computational Physics'' '''390''', 51&nbsp;&ndash;&nbsp;65.
<br/>

<div id="Lehmkuhl2007"></div>
Lehmkuhl, O., Perez Segarra, C.D., Borrell, R., Soria, M. & Oliva, A. 2007 
:Termofluids: A new Parallel unstructured CFD code for the simulation of turbulent industrial problems on low cost PC cluster. ''Proceedings of the Parallel CFD Conference'' pp. 1&nbsp;&ndash;&nbsp;8.
<br/>

<div id="Lilly1992"></div>
Lilly, D. K. 1992 
:A proposed modification of the Germano subgrid-scale closure method. ''Physics of Fluids A'' '''4''' (3), 633&nbsp;&ndash;&nbsp;635.
<br/>

<div id="Lin2007"></div>
Lin, C.-L., Tawhai, M.H., Mclennan, G. & Hoffman, E.A. 2007 
:Characteristics of the turbulent laryngeal jet and its effect on airflow in the human intra-thoracic airways. ''Respiratory Physiology and Neurobiology'' '''157''', 295&nbsp;&ndash;&nbsp;309.
<br/>

<div id="Lizal2015"></div>
Lizal, Frantisek, Belka, Miloslav, Adam, Jan, Jedelsky, Jan & Jicha, Miroslav 2015 
:A method for ''in vitro'' regional aerosol deposition measurement in a model of the human tracheobronchial tree by the positron emission tomography. ''Proceedings of the Institution of Mechanical Engineers, Part H: Journal of Engineering in Medicine'' '''229''' (10), 750&nbsp;&ndash;&nbsp;757.
<br/>

<div id="Lizal2012"></div>
Lizal, Frantisek, Elcner, Jakub, Hopke, Philip K, Jedelsky, Jan & Jicha, Miroslav 2012 
:Development of a realistic human airway model. ''Proceedings of the Institution of Mechanical Engineers, Part H: Journal of Engineering in Medicine'' '''226''' (3), 197&nbsp;&ndash;&nbsp;207.
<br/>

<div id="Matida2004"></div>
Matida, E. A., Finlay, W. H., Lange, C. F. & Grgic, B. 2004 
:Improved numerical simulation of aerosol deposition in na idealized mouth-throat. ''Aerosol Science'' '''35''', 1&nbsp;&ndash;&nbsp;19.
<br/>

<div id="Mead-Hunter2013"></div>
Mead-Hunter, Ryan, King, Andrew J.C., Larcombe, Alexander N. & Mullins, Benjamin J. 2013 
:The influence of moving walls on respiratory aerosol deposition modelling. ''Journal of Aerosol Science'' '''64''', 48&nbsp;&ndash;&nbsp;59.
<br/>

<div id="Muela2019"></div>
Muela, J., Rigoala, J., Oliet, C., Perez-Segarra, C.D. & Oliva, A. 2019 
:Assessment of numerical aspects using LES in particle separation devices. ''In 8th European Conference for Aeronautics and Aerospace Sciences (EUCASS)'', p. 678.
<br/>

<div id="Nicoud1999"></div>
Nicoud, Franck & Ducros, Frederic 1999 
:Subgrid-scale stress modelling based on the square of the velocity gradient tensor. ''Flow, turbulence and Combustion'' '''62''' (3), 183&nbsp;&ndash;&nbsp;200.
<br/>

<div id="OpenFOAM2013a"></div>
OpenFOAM Foundation 2013a 
:OpenFOAM Programmer’s Guide, version 2.2.1 edn. London, UK.
<br/>

<div id="OpenFOAM2013b"></div>
OpenFOAM Foundation 2013b
:OpenFOAM User Guide, version 2.2.1 edn. London, UK.
<br/>

<div id="Radhakrishnan2009"></div>
Radhakrishnan, H. & Kassinos, S. 2009 
:CFD modeling of turbulent flow and particle deposition in human lungs. ''31st Annual International Conference of the IEEE EMBS'', Mineapolis, Minnesota, USA pp. 2867&nbsp;&ndash;&nbsp;2870.
<br/>

<div id="Sagaut2006"></div>
Sagaut, P. 2006 
:Large Eddy Simulation for Incompressible Flows: An Introduction. ''Springer''.
<br/>

<div id="Schmidt2004"></div>
Schmidt, Andreas, Zidowitz, Stephan, Kriete, Andres, Denhard, Thorsten, Krass, Stefan & Peitgen, Heinz-Otto 2004 
:A digital reference model of the human bronchial tree. ''Computerized Medical Imaging and Graphics'' '''28''' (4), 203&nbsp;&ndash;&nbsp;211.
<br/>

<div id="Sciacchitano2016"></div>
Sciacchitano, A. and Wieneke, B. 2016
:PIV uncertainty propagation. ''Meas. Sci. Technol.'' '''27''' (084006), 20pp.
<br/>
<div id="Smagorinsky1963"></div>
Smagorinsky, Joseph 1963 
:General circulation experiments with the primitive equations: I. The basic experiment. ''Monthly weather review'' '''91''' (3), 99&nbsp;&ndash;&nbsp;164.
<br/>

<div id="Tabor2004"></div>
Tabor, G. R., Baba-Ahmadi, M. H., de Villiers, E. & Weller, H. G. 2004 
:Construction of inlet conditions for LES of turbulent channel flow. ''European Congress on Computational Methods in Applied Sciences and Engineering (ECCOMAS)''.
<br/>

<div id="Tawhai2011"></div>
Tawhai, M.H. & Lin, C.-L. 2011 
:Airway gas flow. ''Comprehensive Physiology'' '''1''', 1135&nbsp;&ndash;&nbsp;1157. 
<br/>

<div id="Trias2011"></div>
Trias, F. X. & Lehmkuhl, O. 2011
:A self-adaptive strategy for the time integration of Navier-Stokes equations. ''Numerical Heat Transfer. Part B'' '''60''' (2), 116&nbsp;&ndash;&nbsp;134.
<br/>
<div id="Vazquez2016"></div>
Vazquez, A. M., Houzeaux, G., Koric, S., Artigues, A., Aguado-Sierra, J., Aris, R., Mira, D., Calmet, H., Cucchietti, F., Owen, H., Casoni, E., Taha, A., Burness, E. D., Cela, J. M. & Valero, M. 2016 
:Alya: Multiphysics engineering simulation towards exascale. ''Journal Computational Sciences'' '''14''', 15&nbsp;&ndash;&nbsp;27.
<br/>

<div id="Verstappen2011"></div>
Verstappen, Roel 2011 
:When does eddy viscosity damp subfilter scales sufficiently? ''Journal of Scientific Computing'' '''49''' (1), 94.
<br/>

<div id="Verstappen2010"></div>
Verstappen, R.W.C.P., Bose, S.T., Lee, J., Choi, H. & Moin, P.P 2010 
:A dynamic eddy-viscosity model based on the invariants of the rate-of-strain. ''In Proceedings of the summer program'', pp. 183&nbsp;&ndash;&nbsp;192. Center for Turbulence Research, Stanford University Stanford.
<br/>

<div id="Verstappen2003"></div>
Verstappen, R.W.C.P. & Veldman, A.E.P. 2003 
:Symmetry-preserving discretization of turbulent flow. ''Journal of Computational Physics'' '''187''', 343&nbsp;&ndash;&nbsp;368.
<br/>

<div id="Weibel1963"></div>
Weibel, E. R. 1963 
:Morphometry of the human lung. Springer-Verlag, Berlin.
<br/>

<div id="Wu2014"></div>
Wu, Dan, Tawhai, Merryn H, Hoffman, Eric A & Lin, Ching-Long 2014 
:A numerical study of heat and water vapor transfer in MDCT-based human airway models. ''Ann Biomed Eng'' '''42''' (10), 2117&nbsp;&ndash;&nbsp;2131.
<br/>

<div id="Xi2018"></div>
Xi, Jinxiang, April Si, Xiuhua, Dong, Haibo & Zhong, Hualiang 2018 
:Effects of glottis motion on airflow and energy expenditure in a human upper airway model. ''European Journal of Mechanics - B/Fluids'' '''72''', 23&nbsp;&ndash;&nbsp;37.
<br/>

<div id="Yin2010"></div>
Yin, Youbing, Choi, Jiwoong, Hoffman, Eric A., Tawhai, Merryn H. & Lin, Ching-Long 2010 
:Simulation of pulmonary air flow with a subject-specific boundary condition. ''Journal of Biomechanics'' '''43''' (11), 2159&nbsp;&ndash;&nbsp;2163.
<br/>

<div id="Zang1993"></div>
Zang, Yan, Street, Robert L. & Koseff, Jeffrey R. 1993 
:A dynamic mixed subgrid-scale model and its application to turbulent recirculating flows. ''Physics of Fluids A: Fluid Dynamics'' '''5''' (12), 3186&nbsp;&ndash;&nbsp;3196.
<br/>

<div id="Zhang2004"></div>
Zhang, Z. & Kleinstreuer, C. 2004 
:Airflow structures and nano-particle deposition in a human upper airway model. ''Journal of Computational Physics'' '''198''', 178&nbsp;&ndash;&nbsp;210.
<br/>

<div id="Zhou2005"></div>
Zhou, Yue & Cheng, Yung-Sung 2005 
:Particle deposition in a cast of human tracheobronchial airways. ''Aerosol Science and Technology'' '''39''' (6), 492&nbsp;&ndash;&nbsp;500.


<br/>
----

{{ACContribs
|authors=P. Koullapis<sup>a</sup>, J. Muela<sup>b</sup>, O. Lehmkuhl<sup>c</sup>, F. Lizal<sup>d</sup>, J. Jedelsky<sup>d</sup>, M. Jicha<sup>d</sup>, T. Janke<sup>e</sup>, K. Bauer<sup>e</sup>, M. Sommerfeld<sup>f</sup>, S. C. Kassinos<sup>a</sup>
|organisation=<br><sup>a</sup>Department of Mechanical and Manufacturing Engineering, University of Cyprus, Nicosia, Cyprus<br>
<sup>b</sup>Heat and Mass Transfer Technological Centre, Universitat Politècnica de Catalunya, Terrassa, Spain<br>
<sup>c</sup>Barcelona Supercomputing center, Barcelona, Spain <br>
<sup>d</sup>Faculty of Mechanical Engineering, Brno University of Technology, Brno, Czech Republic<br>
<sup>e</sup>Institute of Mechanics and Fluid Dynamics, TU Bergakademie Freiberg, Freiberg, Germany<br>
<sup>f</sup>Institute Process Engineering, Otto von Guericke University, Halle (Saale), Germany <br>
}}
{{ACHeader
|area=7
|number=02
}}

© copyright ERCOFTAC 2020