氣溶膠粒子在呼吸系統(tǒng)內(nèi)沉降規(guī)律的數(shù)值模擬

趙 卓，黃 平

氣溶膠粒子在呼吸系統(tǒng)內(nèi)沉降規(guī)律的數(shù)值模擬

趙 卓，黃 平

目的：研究可吸入氣溶膠粒子在人體氣管、支氣管內(nèi)的沉積規(guī)律。方法：建立基于Weibel模型的G0～G3級氣管、支氣管模型，在呼吸強度為15、30和60 L/min時，設定不同的出口邊界條件，模擬氣溶膠粒子不同密度和不同粒子直徑時的顆粒沉積情況。結果：不同的邊界條件下，流場變化不大，對粒子的沉積影響不大。相對于粒子密度，粒子的直徑和呼吸強度對沉積規(guī)律的影響更大。結論：該研究對于開發(fā)新的氣溶膠吸入治療藥物和新的氣溶膠吸入治療裝置具有一定的借鑒意義，為進一步研究氣溶膠濃度和失能劑起效的血藥濃度的對應關系提供了一定的依據(jù)。

氣溶膠粒子；呼吸系統(tǒng)；氣管支氣管；沉積

0 引言

隨著局部武裝沖突與全球恐怖活動的不斷增加，傳統(tǒng)意義上的彈藥已不能滿足日益復雜的安全局面。使人暫時失去抵抗能力而不致死的非致命武器應運而生，其中，化學類非致命武器是利用機械噴灑、燃燒、爆炸等氣霧釋放裝置將催淚劑、麻醉劑、失能劑等化學物質(zhì)分散成可吸入的氣溶膠粒子，通過人的呼吸系統(tǒng)進入人體，其起效時間短且持續(xù)時間長，使人迅速喪失戰(zhàn)斗力，最終達到控制人員的目的。因此，氣溶膠粒子在人體呼吸系統(tǒng)內(nèi)沉積規(guī)律的研究成為一大熱點[1]。

氣溶膠粒子在人體呼吸系統(tǒng)內(nèi)的沉積規(guī)律的研究方法目前主要有2種：實物實驗和數(shù)值模擬。其中，實物實驗成本較高、可重復性差，且操作困難。而借助計算機的流體動力學數(shù)值模擬（computational fluid dynamics simulation，CFD）經(jīng)濟實用、可重復性好。越來越多的研究者開始運用數(shù)值模擬方法來研究可吸入氣溶膠粒子的沉積規(guī)律。運用CFD方法，對流場計算的湍流模型的選取有標準k-ε模型、重整化群（renormalization group，RNG）k-ε模型、標準k-ω模型、大渦模擬（large eddy simulation，LES）等。Staple－ton等[2]采用標準k-ε模型模擬了多分散的氣溶膠粒子在理想的口喉模型中的沉積情況，并和實驗進行了對比，指出標準k-ε模型在模擬氣溶膠粒子沉積方面并不是很準確。曾敏捷等[3]采用RNG k-ε模型模擬了氣溶膠粒子在完整上呼吸道內(nèi)的沉積情況，與Grgic[4]的實驗結果對比表明，RNG k-ε模型可以較好地模擬管流，模擬結果比較準確。Zhang等[5]采用標準k-ω模型研究了微納米粒子在人體上呼吸道內(nèi)的沉積情況，證實可靠的數(shù)值模擬方法可為氣溶膠吸入治療提供指導。Luo等[6]采用LES對非對稱的氣管支氣管模型進行了數(shù)值模擬，結果表明LES方法是可靠且可行的，如果氣管支氣管模擬能更接近人體真實的氣管支氣管，LES方法則可以更好地預測氣溶膠粒子的沉積。

對人體呼吸系統(tǒng)內(nèi)可吸入氣溶膠粒子沉積的研究，邊界條件的選取非常重要。入口邊界選用第一類邊界條件，即速度入口邊界，而出口邊界條件的選取有爭議。在目前的研究中，多數(shù)采用壓力出口（pressure-outlet），還有一部分采用出流邊界（outflow）類型。周鑫[7]認為采用自由出流的邊界類型誤差比較大。Nowak[8]在研究Weibel模型時，對比了壓力出口和自由出流邊界2種類型，從氣流流場和粒子沉積率2個方面進行了數(shù)值模擬，結果沒有發(fā)現(xiàn)顯著的差異。本文運用CFD方法，選取標準k-ω模型，采用壓力出口和自由出流出口2種邊界類型，研究了可吸入氣溶膠粒子在人體支氣管樹中的G0～G3內(nèi)的沉積情況。

1 模型描述

1.1 人體氣管支氣管幾何模型

本文采用廣泛使用的Weibel[9]模型，氣管支氣管簡化為一系列光滑的圓形管道，同一等級的支氣管是對稱的。詳細數(shù)據(jù)見表1，氣管支氣管模型圖如圖1所示。

表1 G0～G3氣管支氣管幾何參數(shù)

圖1 氣管支氣管模型圖

1.2 控制方程

根據(jù)人體呼吸道的結構特點，本文采用標準的k-ω模型[10]對流場進行計算，其輸運方程為

連續(xù)性方程：

動量方程：

湍動能（k）和耗散率（ω）由以下2個方程解得：

式（1）～（5）中，ui和xj（i和j取1、2、3）分別表示x、y、z方向上的速度分量和空間坐標分量，t、ρ、p、ν、τij分別表示時間、密度、壓力、運動黏度和雷諾應力張量，其單位分別為s、kg/m3、N、m2/s、Pa，而cμ、α、β、β*、σk、σω則為湍流常量。

而對于顆粒的運動，本文選用基于牛頓第二定律的拉格朗日離散相模型，通過分析氣溶膠顆粒的受力情況，從而得到顆粒的軌跡。顆粒的受力可能包括布朗力、Basset力、Magnus升力、Staffman升力、熱泳力、電泳力、重力、曳力等。本文研究的粒子在1～10 μm之間，所以，只考慮重力、曳力和Staffman升力，其中，曳力是最主要的。顆粒的受力方程為

式中，右邊3項分別為顆粒受到的Stokes阻力、重力和Staffman升力，upi為顆粒在i方向的速度分量，CD為阻力系數(shù)。

1.3 邊界條件和數(shù)值計算方法

（1）入口邊界：入口氣流為常溫狀態(tài)下的不可壓縮的空氣，密度ρ=1.025 kg/m3，黏度μ=18.1×10-6Pa·s，采用速度入口邊界（velocity-inlet），給定穩(wěn)態(tài)吸氣時入口截面方向的平均速度值。

（2）出口邊界：本文選用2種出口邊界條件，分別是壓力出口和出流邊界條件。

（3）固壁邊界：氣管支氣管的壁面采用無滑移的固壁邊界，也不考慮壁面與流體的摩擦而產(chǎn)生的熱效應。加入顆粒相后固壁邊界為“trap”，即捕獲。

對于氣管支氣管模型的離散化，對整個計算區(qū)域分別建立了40萬、70萬、100萬的網(wǎng)格。經(jīng)過大量的數(shù)值模擬計算，70萬的網(wǎng)格數(shù)量可以滿足精度的需要，計算所耗資源合理，因此，本文采用約70萬的網(wǎng)格單元，部分網(wǎng)格如圖2所示。本文采用Fluent商業(yè)軟件，選取分離隱式求解納維—斯托克斯（Navier-Stokes）方程組，壓力速度耦合為SIMPLEC，壓力、動量、湍流動能和耗散率的離散格式采用二階迎風格式。

1.4 模型驗證

圖2 局部網(wǎng)格

為了對本文模型和計算方法的正確性進行驗證，選取15和30L/m in進行模擬，并與實驗結果[11-12]進行對比，結果如圖3所示。從圖3中的數(shù)據(jù)可知，數(shù)值計算結果和原有實驗數(shù)據(jù)能夠很好地吻合，誤差在允許的范圍內(nèi)，表明本文所建立的氣管支氣管模型和所采用的數(shù)值仿真方法是合理的，可以用于可吸入氣溶膠粒子沉積規(guī)律的研究。

圖3 實驗與數(shù)值模擬的對比

2 結果和分析

2.1 氣相流場

從圖4和圖5中可以看出，不同的呼吸流量氣管支氣管內(nèi)的速度云圖是相似的，呼吸流量的增加沒有影響氣流的結構，只是增大了支氣管下游的速度最大值。最大速度值和呼吸流量呈現(xiàn)正相關性，采用壓力出口的流量為15、30、6 0 L/m in時，速度最大值分別為1.9、3.4、6.5 m/s。而采用出流出口邊界的速度最大值分別為 2、3.8、7 m/s，均稍高于對應的壓力出口邊界的速度最大值。從圖中可以看出，支氣管中的氣流比較平穩(wěn)，但在上下兩級氣管的交接處，由于受到較大的阻力，氣流變化較大。在同一級支氣管中，靠近上級支氣管中心軸線處的速度較大。

圖4 壓力出口

2.2 氣溶膠顆粒沉積的影響因素分析

2.2.1 呼吸流量的影響

從圖6可以看出，呼吸流量是影響粒子沉積的一個重要因素。1 μm的粒子，呼吸流量從15 L/m in增加到7 5 L/m in，壓力出口和出流出口的粒子沉積率分別增加32.52%和32%；5 μm的粒子，呼吸流量從15 L/m in增加到7 5 L/m in，壓力出口和出流出口的粒子沉積率分別增加32.1%和28.1%；9 μm的粒子，呼吸流量從15 L/m in增加到7 5 L/m in，壓力出口和出流出口的粒子沉積率分別增加 37.3%和37.1%。從圖上還可以看出，對于1 μm的粒子，呼吸流量由15 L/m in變?yōu)?0 L/m in時，粒子的沉積率幾乎沒有發(fā)生變化，表明對于粒徑較小的粒子，呼吸流量較低時，增加呼吸流量，對粒子的沉積影響不大。對于1、3和5 μm的粒子，呼吸流量從15 L/m in增加到30 L/m in時，粒子的沉積率變化不大。但當呼吸流量上升為7 5L/m in時，粒子的沉積率變化明顯。

圖5 出流出口

2.2.2 粒徑的影響

從圖7可以看出，粒徑對可吸入氣溶膠粒子的影響很明顯。粒徑較小的顆粒因慣性較小而跟隨流體運動，不易吸附在氣管支氣管壁面上；而隨著粒徑的增大，慣性隨之增大，更容易附著在壁面上。呼吸流量為15 L/m in時，粒徑從1～9 μm，壓力出口和出流出口的粒子沉積率分別增加 27.02%和28.30%；呼吸流量為30 L/m in時，粒徑從1～9 μm，壓力出口和出流出口的粒子沉積率分別增加45.19%和46.96%；呼吸流量為6 0 L/m in時，粒徑從1～9 μm，壓力出口和出流出口的粒子沉積率分別增加41.9%和36.5%。

圖6 沉積率和流量的關系

圖7 沉積率和粒徑的關系

2.2.3 密度的影響

如圖8所示，本文選取了Q=15 L/min時，粒徑為1、5、9 μ m的粒子在密度為1 000、1 250、1 500、1 750、2 000 kg/m3時的粒子沉積情況?？傮w來看，粒子的沉積率隨著密度的增加而增加。出口邊界不同時，粒子的沉積率隨密度變化不大。當粒子直徑為1 μ m時，氣溶膠粒子的沉積率隨密度的變化不大；但是當粒徑為9 μ m時，氣溶膠粒子的沉積率增加較大。對于粒徑大的粒子，粒子的沉積機制主要由慣性力主導，而慣性力和密度呈正相關。因此，對于氣溶膠吸入治療來說，不同粒徑的粒子在不同的呼吸強度下，密度不同時沉積率不同。

圖8 沉積率和密度的關系

3 討論

本文研究了穩(wěn)態(tài)吸氣下可吸入氣溶膠粒子的人體氣管支氣管內(nèi)的運動沉積規(guī)律，分析了氣管支氣管內(nèi)的氣流運動形式和變化規(guī)律，對比了壓力出口和出流出口2種邊界條件類型，得到以下結論：

（1）壓力出口和出流出口邊界都可以合理地反映人體呼吸道內(nèi)的流場狀態(tài)，但對于呼吸道內(nèi)的氣場，這2種邊界類型的影響不大。對于氣管支氣管級數(shù)比較少的，采用壓力出口和出流出口邊界對流場和氣溶膠粒子的沉積影響不大。

（2）影響可吸入氣溶膠粒子的因素有粒子的密度、動力學直徑、呼吸流量等，其中，粒子的動力學直徑和呼吸流量是有較大影響的2個因素。呼吸強度和氣溶膠粒子的慣性是影響粒子沉積的主要因素。在研究可吸入藥物時，對確定的藥物粒子，密度是個定值，應該根據(jù)所要達到的效果和想要粒子沉積的位置，確定一個合理的粒子直徑和進氣速度。

本文的研究對于開發(fā)新的氣溶膠吸入治療藥物和新的氣溶膠吸入治療裝置具有一定的借鑒意義，可為進一步研究氣溶膠濃度和失能劑起效的血藥濃度的對應關系提供一定的依據(jù)。同時，還有很多工作要完善，比如建立完整的包括口、咽、喉、氣管和多級支氣管的人體呼吸道模型；建立兒童的呼吸道模型，研究氣溶膠粒子在兒童呼吸道內(nèi)的沉積情況；考慮吸入裝置的實際模型，細化口腔入口處的速度和粒子入射情形等。

[1]Kieinstreuer C.A new methodology for targeting drug-aerosols in the human respiratory system[J].International Journal of Heat and Mass Transfer，2008，51（23-24）：5 578-5 589.

[2]Stapleton K W.On the suitability of turbulence modeling for aerosol deposition in the mouth and throat：a comparison with experiment[J]. Journal of Aerosol Science，2000，31（6）：739-749.

[3]曾敏捷，胡桂林，樊建人.微顆粒在人體上呼吸道中運動沉積的數(shù)值模擬[J].浙江大學學報：工學版，2006，40（7）：1 164-1 167，1 256.

[4]Grgic B，F(xiàn)inlay W H，Heenan A F.Regional aerosol deposition and flow measurements in an idealized mouth and throat[J].Journal of Aerosol Science，2004，35（1）：21-32.

[5]ZHANG Zhe，Kleinstreuer C，Donohue J F，et al.Comparison of micro-and nano-size particle depositions ina human upper airway model[J].Journal of Aerosol Science，2005，36（2）：211-233.

[6]LUO X Y.LES modelling of flow in a simple airway model[J].Medical Engineering&Physics，2004，26（5）：403-413.

（????）（????）

[7]周鑫.顆粒物在人體呼吸系統(tǒng)中傳輸與沉積的數(shù)值模擬研究[D].長沙：中南大學，2010.

[8]Nowak N，Kakade P P，Annapragada A V.Computational fluid dynamics simulation of airflow and aerosol deposition in human lungs[J]. Annals of Biomedical Engineering，2003，31（4）：374-390.

[9]Weibel E R.Morphometry of the human lung[M].New York：Springer，1963：1-3.

[10]LI Z，Kieinstreuer C，ZHANG Z.Particle deposition in the human tracheobronchial airways due to transient inspiratory flow patterns[J]. Journal of Aerosol Science，2007，38（6）：625-644.

[11]ZHANG Z，Kieinstreuer C，Kim C S.Comparison of analytical and CFD models with regard to micron particle deposition in a human 16-generation tracheobronchial airway model[J].Journal of Aerosol Science，2009，40（1）：16-28.

[12]Emmett P，Aitken R，Hannan W.Measurements of the total and regional deposition of inhaled particles in the human respiratory tract[J]. Journal of Aerosol Science，1982，13（6）：549-560.

（收稿：2013-10-24 修回：2014-06-23）

Numerical simulation on deposition of aerosol particles in respirator trace

ZHAO Zhuo,HUANG Ping
(State Key Laboratory of Explosion Science and Technology,Beijing Institute of Technology,Beijing 100081,China)

ObjectiveTo study the deposition of inhalable aerosol particles in the trachea and bronchia.MethodsA tracheal bronchus model of generation G0 to G3 based on Weibel's model was built.The aerosol particles deposition was simulated under breathing intensity of 15,30 and 60 L/min,with different outlet boundary conditions,densities and diameters.ResultsFlow field differed little as the boundary changed,while the breath densities and particle diameters had more influences on the deposition.ConclusionThe study can help to develop aerosol inhalable medicine and apparatus, and to determine the correlation between the aerosol concentration and blood drug level when the incapacitating agent comes into effect.[Chinese Medical Equipment Journal，2014，35（10）：21-24，31]

aerosol particle;respiration system;tracheal bronchus;deposition

R318；TU196.2

A

1003-8868（2014）10-0021-05

10.7687/J.ISSN1003-8868.2014.10.021

趙 卓（1987—），男，研究方向為功能材料，E-mail：2120110343@ bit.edu.cn。

100081北京，北京理工大學爆炸科學與技術國家重點實驗室（趙 卓，黃 平）