WORKGROUP FOR

MULTIPHAS FLOWS

Ansprechpartner:

Prof. Dr.-Ing. Habil Martin Sommerfeld

1. Introduction

The understanding of fine particle transport and deposition in the upper respiratory tract of the lung is important for delivery of aerosolised drug powders as well as an unwanted exposure to environmental aerosols. The breathing flow within the lung is rather complex and of course strongly time-dependent. In addition, the tissue is not rigid and has some flexibility which is also varying in time. An experimental analysis of fine particle deposition especially in vitro is reasonably sophisticated, although conducted with a number of simplifications (Lizal et al. 2015). Therefore, computational fluid dynamics has been developed over the past decades for being an efficient analysis tool for fine particle deposition in lungs. Naturally, an Euler/Lagrange approach (Sommerfeld et al. 2008, Sommerfeld 2017) is the best choice for doing such computations (for applications see for example Longest and Xi 2007, Kolanjiyil and Kleinstreuer 2016). Mostly the computations assume stationary flow at pre-defined typical flow rates. There are of course several options for simulating the flow through lung models with regard to turbulence modelling (Koullapis et al. 2018). Large eddy simulations (LES), using rather simple sub-grid scale turbulence models, are probably the best choice but rather expensive as relatively fine grids are needed and computations of flow and particles have to be done simultaneously and unsteady. On the other hand, RANS methods (Reynolds-averaged Navier-Stokes) are numerically much faster, but require more or less sophisticated turbulence models for predicting such flows with close to laminar and turbulent flow regions and associated transitions.

A key requirement for analysing different turbulence models for predicting respiratory systems is of course a validation based on experimental data. Detailed PIV (particle image velocimetry) measurements in a rather complex lung model (comprising the oral cavity, larynx and tracheobronchial airways down to the 7th generation of branching) using a refractive index matching technique were recently conducted by Janke et al. (2019) providing mean velocity and turbulent kinetic energy fields and profiles in different regions of the lung model. This test case was the basis for analysing the performance of different RANS turbulence models and the influence of inlet conditions. Based on the RANS flow field Lagrangian particle tracking was applied for analysing the consequences of the predicted turbulence field on the regional and size-dependent deposition of the fine particles. Local particle deposition measurements were conducted in the same model and were presented by Lizal et al. (2015). Both experimental cases namely flow and deposition measurements are included in the ERCOFTAC data base for further studies (see Koullapis et al. 2020 and Koullapis et al. 2019).

2. RANS Results

The present computations were done using a RANS-based Euler/Lagrange approach, however without two-way coupling, on the basis of the OpenFOAM package which was thoroughly extended to allow proper particle tracking computations (Sgrott and Sommerfeld 2019). Several turbulence models were considered, namely, the k-w-SST, standard k-e type models and the Reynolds stress (RSM-SSG) turbulence model. The geometry of the lung model is illustrated in Figure 1 a) with inflow (diameter 20 mm), oral cavity, trachea and the bronchia region. A short straight inlet pipe (here 30 mm length) is considered before the oral cavity and the 10 outlets with zero-gradient conditions are indicated by arrows (Koullapis et al. 2018). After preliminary studies the lung model was discretise with a grid of about 7.9 million cells having a refined region near the wall for all the presented RANS computations. Besides the analysis of different turbulence models also the inlet boundary condition was varied considering the following cases:

• Mapped inlet over a distance of 20 mm of the short inlet pipe shown in Figure 1 a)
• Plug flow inlet condition at the entrance of the short inlet pipe with a turbulence intensity of 5%
• A 10-diameter extension of the short inlet pipe (see Figure 1 a) with plug flow inlet and a 5% turbulence intensity.

For the outlet all variables (except for the pressure) were set to zero gradient, and at the walls the different wall functions were applied (i.e. for turbulent kinetic energy, Reynolds-stresses, dissipation rate, specific dissipation rate and turbulent viscosity). For the pressure a fixed value of zero was specified at the outlet and zero gradient was used at the inlet.

From two computations vertical cuts through the upper part of the lung model are presented in Figure 1 b). Both simulations show roughly the same mean flow structure, however the k-e flow field presents less fine flow structures which are better reflected in the RSM-SSG results and low velocity regions show steeper edges and gradients. The flow field obtained with the k-w-SST turbulence model is very similar to the k-e result however provides a bit stronger shear layers. It should be emphasised RSM simulations were not completely stable and not fully stationary even after several thousands of iterations.

Very large differences are however obtained for the predicted turbulent kinetic energy fields when using the three turbulence models considered (see Figure 2). The lowest levels are obtained with the k-w-SST and the Reynolds-stress model (RSM-SSG) yields expectedly the highest values which are also in better agreement with the measurements being discussed below. Surprisingly the standard k-e turbulence model yields higher turbulence levels as the k-w-SST model.

3. Particle Deposition

The dispersed phase is treated in the Lagrangian framework, in which the particles are simultaneously and time-dependent tracked through the stationary flow. The particles are treated as point-masses and assumed to be spherical. Only one-way coupling simulations are considered, since the particle volume fraction is very small to have a relevant effect on the continuous phase. The positions and velocities of the particles are calculated by integrating the equations of motion considering the drag force with Cunningham correction, gravity/buoyancy, the slip-shear lift force and the random Brownian force (Sommerfeld et al. 2008). Since the Brownian force is under the present conditions only relevant for particle smaller than 0.5 mm it was omitted for the particle sizes considered here. Moreover, drag and lift forces were extended to account for wall effects (Zeng et al. 2009) which was shown to slightly reduce deposition rate (Sgrott and Sommerfeld 2019). Since the considered particles are normally below 10 mm any rotational effects can be neglected. All particles hitting the wall are assumed to be deposited. The time step of the particle tracking calculation is automatically and independently adjusted along the trajectories by considering all relevant time scales which are also changing throughout the flow field, namely the particle response time  and the integral time scale of turbulence  obtained from the turbulence model.

In order to model the effect of turbulence on particle transport and eventually deposition, the instantaneous fluid velocity has to be used in the fluid forces. From the calculation of the fluid field only the mean velocities are available, wherefore the instantaneous fluid fluctuating velocities have to be generated with the help of the calculated turbulence properties. For that purpose, the single-step Langevin equation model is used (Sommerfeld et al. 1993, Sommerfeld et al. 2008). Unfortunately these kind of turbulent dispersion models are subject to a spurious drift (Sommerfeld et al. 1993, Sommerfeld 1996) which will cause an unrealistic enhancement of fine particle deposition (Matida et al. 2003). There are several approaches for correcting such a spurious drift (Legg and Raupach 1982, Dehbi 2008), but for two-equation turbulence models they are approximately. Only if the full Reynolds stress tensor is considered for the drift correction with an RSM model this effect may be completely corrected.

4. Results and discussion

Numerical computations were conducted for different inlet and outlet conditions, a number of flow rates and various particles sizes. The computed profiles of flow velocity and turbulent kinetic energy are compared with measurements (Janke et al. 2019) and also with high resolution LES using 50 million grids (Koullapis et al. 2019). First the influence of the different inlet conditions is considered. For that purpose, a vertical profile in the oral cavity (Profile B in Figure 1) was selected which is located just downstream of the inlet. Hence, the profiles of the normalised velocity magnitude and the turbulent kinetic energy are quite sensitive with respect to the inlet conditions, especially at this location close to the inlet. The profile of the velocity magnitude is characterised by a moderate shear layer resulting from the inlet jet and a low velocity region in the upper part of the cavity (see also Figure 1 b). The velocity magnitude is captured reasonably well for all inlet conditions (Figure 3 a). Only the k-w-SST turbulence model with mapped inlet yields a distinct minimum, which is not found in the experiments as well as the LES results. Larger differences are observed in the profiles of the turbulent kinetic energy. Slightly larger values are obtained when using the standard k-e turbulence model with both mapped inlet and specifying a plug flow at the short inlet. The highest turbulence levels are predicted with the k-e model and a 10-diameter long inlet pipe. The resulting values are even quite close to the LES-results at this location. However, such an inlet condition is completely different from the experimental condition using a short inlet pipe (Figure 1) and sucking the liquid out of a surrounding container through the entire lung model (Janke et al. 2019).

Profile located downstream of the oral cavity in the middle region of the trachea (profile D in Figure 1) are compared in Figure 4 using different turbulence models with the same inlet condition, namely, a short inlet pipe with mapped fluid properties over a distance of 20 mm. The mean velocity profile resulting from the RSM computation is in quite good agreement with measurements and LES-results (Figure 4). However, both the k-w-SST and the k-e turbulence models do not correctly capture the location of the shear layer and under-predict the turbulent kinetic energy compared to a high-resolution LES and the measurements. The k-values obtained with the RSM model yields a bit higher values, but still only about 50% from the LES results. However, the shape of the k-profile obtained with the RSM model is very similar to the LES. Naturally such an under-prediction of turbulence is not a good basis for doing fine particle deposition studies since turbulence is a key mechanism for particle deposition.

Particle deposition was analysed along the entire lung model for obtaining the local deposition fraction as well as the total deposition in dependence of particle size (Figure 5). Lagrangian computations were conducted for different force formulations (here only wall corrections are used for drag and lift) and also with and without turbulent particle dispersion model (TD) as well as with full drift correction for the RSM turbulence model. The kind of inlet boundary condition did not have a pronounced influence on the overall deposition and therefore only the short inlet with mapped inlet boundary for the gas flow is considered. The spread of the results is quite large and the LES-values actually have the lowest deposition for the fine particles and the biggest one for the largest particles, i.e. 10 mm. Here it should be emphasised that these simulations do not account for the (SGS) sub-grid-scale turbulent dispersion. For the three results with the different turbulence models including wall effects and standard turbulent dispersion (no drift correction) the magnitudes of particle deposition depend on the intensity of the predicted turbulence, which grows in order for the k-w-SST, the k-e and the RSM-SSG turbulence models. Actually the result for the k-w-SST model is very close to the deposition obtained with the RSM but without turbulent particle transport, even a bit lower. Note that turbulence expectedly enhances deposition. With the drift correction the RSM provides slightly higher deposition fractions for the larger particles compared to the LES. On the other hand for the smaller particles the deposition fraction with the RSM is predicted to be larger than the LES result, which might be an effect of the neglected SGS transport of such fine particles. The available measurement point for 4.3 mm particles (Lizal et al. 2015, Koullapis et al. 2019) is just underneath the LES and just on the line of the RSM-deposition without turbulent transport.

For the Reynolds-stress turbulence model (RSM) the local deposition along the lung model is shown in Figure 6 and compared with measurements and also the LES results. The LES result is in reasonable agreement with the measurements, i.e. the trends are predicted correctly. It should be noted that the numbers of the horizontal axis are not in sequence with respect to the flow passing through the lung model. They rather follow branches and then jump to the other side (see Koullapis et al. 2019). Typical elements are: 3, 5 and 8 are branches, 4 is a kind of curved pipe and 9 as well as 12 are again branches. Such branches should give a bit higher deposition since particle inertia effects become important, e.g. high deposition at 6 and 12. The numbers from 13 to 22 refer to the multiple branches just in front of the outlet buckets. They all have relatively high deposition fraction.

Especially in the initial parts (i.e. from 1 to 12) the RSM with full drift correction shows quite high deposition fractions compared to the measurements and the LES. From 5 to 12 the values of LES are basically constant just as the full RSM, which are however higher. The deposition for the RSM without particle turbulent dispersion follows in the initial part the trends of the LES, being however just a bit lower. Hence, for such a complex geometry it is very hard to accurately predict the correct local deposition since there are so many things which could deviate between experiment and numerical predictions, such as differences in mean flow field as well as turbulence. Very interesting is however, the deposition minimum between 13 and 14, as well as the maximum at around 19, are captured by the simulations.

5. Acknowledgements

The research work was conducted though cooperation within the 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. The authors also would like to thank Dr. Pantelis Koullapis from University of Cyprus for providing the high-resolution LES data and for helpful discussions with respect to setting up the simulation cases. Moreover, we thank Dr. Katrin Bauer from TU Bergakademie Freiberg for providing the experimental data.

References

Dehbi, A.: Turbulent particle dispersion in arbitrary wall-bounded geometries: A coupled CFD-Langevin-equation based approach. International Journal of Multiphase Flow, Volume 34, 819 – 828, (2008)

Janke, T., Koullapis, P., Kassinos, S.C. and Bauer, K.: PIV measurements of the SimInhale benchmark case. European J. Pharmaceutical Sciences, Vol. 133, 183 – 189 (2019)

Kolanjiyil, A.V. and Kleinstreuer, C.: Computationally efficient analysis of particle transport and deposition in a human whole-lung-airway model. Part I: Theory and model validation. Computers in Biology and Medicine, Vol. 79, 193 – 204 (2016)

Koullapis, P., Kassinos, S.C., Muela, J., Segarra, C., Rigola, J., Lehmkuhl, O., Cui, Y., Sommerfeld, M., Elcner, J., Jisha, M., Saveljic, I., Filipovic, N. Lizal, F. and Nicolaou, L.: Regional aerosol deposition in the human airways: The SimInhale benchmark case and a critical assessment of in silico methods. European Journal of Pharmaceutical Sciences, Vol. 113, 77 – 94 (2018)

Koullapis, P., Lizal, F., Jedelsky, J., Nicolaou, L., Sgrott, O., Jicha, M., Sommerfeld, M., and Kassinos, S.C.: Aerosol deposition in the human upper airways. ERCOFTAC App. Challenge, October 2019

Koullapis, P., Muela, J., Lehmkuhl, O., Lizal, F., Jedelsky, J., Jicha, M., Jahnke, Th., Bauer, K., Sommerfeld, M. and Kassinos, S.C.: Airflow in the human upper airways. ERCOFTAC App. Challenge, July 2020

Lizal, F., Belka, M., Adam, J., Jedelsky, J., Jicha, M.: A method for in vitro regional aerosol deposition measurement in a model of the human tracheobronchial tree by the positron emission tomography. Proc. Inst. Mech. Eng. H J. Eng. Med. Vol. 229, 750–757 (2015)

Legg, B.J. and Raupach, M.R.: Markov-chain simulations of particle dispersion in inhomogeneous flows: The mean drift velocity induced by a gradient in Eulerian velocity variance. Boundary-Layer Meteorology, Vol. 24, 3-13 (1982)

Longest, P.W. and Xi, J.: Effectiveness of direct Lagrangian tracking models for simulating nanoparticle deposition in the upper airways. Aerosol Science and Technology, Vol. 41, 380 – 397 (2007)

Matida, E.A., DeHaan, W.H., Finlay, W. H. and Lange, C.F.: Simulation of Particle Deposition in an Idealized Mouth with Different Small Diameter Inlets. Aerosol Science and Technology, Volume 37, 924-932 (2003)

Sgrott Junior, O.L. and Sommerfeld, M.: Influence of inter-particle collisions and agglomeration on cyclone performance and collection efficiency. Can. J. Chem. Eng. Vol. 97, 511 – 522 (2019)

Sgrott Jr, O.L. and Sommerfeld, M.: Modelling and Simulation of Fine Powder Deposition in Lung Models and Basic Configurations. 10^th International Conference on Multiphase Flow, ICMF 2019, Rio de Janeiro, Brazil, May 19 – 24, 2019

Sommerfeld, M., Kohnen, G. and Rüger, M.: Some open questions and inconsistencies of Lagrangian Particle dispersion models. Ninth Symposium on Turbulent Shear Flows, Kyoto, Paper 15.1. (1993)

Sommerfeld, M.: Modellierung und numerische Berechnung von partikelbeladenen turbulenten Strömungen mit Hilfe des Euler/Lagrange-Verfahrens. Habilitationsschrift, Universität Erlangen-Nürnberg, Shaker Verlag, Aachen (1996)

Sommerfeld, M., van Wachem, B. and Oliemans, R.: Best Practice Guidelines for Computational Fluid Dynamics of Dispersed Multiphase Flows. ERCOFTAC (European Research Community on Flow, Turbulence and Combustion, ISBN 978-91-633-3564-8 (2008).

Sommerfeld, M.: Numerical methods for dispersed multiphase flows. In: Particles in Flows (Eds. T. Bodnár, G.P. Galdi, Š. Necčasová), Series Advances in Mathematical Fluid Mechanics, Springer International Publishing, 327 – 396 (2017)

Zeng, L., Najjar, F., Balachandar, S. and Fischer, P.: Forces on a finite-sized particle located close to a wall in a linear shear flow. Physics of Fluids, Vol. 21, 033302 (2009)