風速對自然通風教室內咳嗽飛沫液滴擴散的影響

曾凱磊 陳燕 譚堰琴 朱清源

東華大學環境科學與工程學院

室內環境中的呼吸道傳染病一直是困擾人們的公共衛生問題，比較典型的呼吸道傳染病有——從發現之日起一直與人類共存的流感和肺結核，曾造成巨大危害但逐步被消滅的SARS 以及2019 年年末爆發的至今仍在世界各地蔓延肆虐的新型冠狀肺炎疫情。呼吸道傳染病的病原體主要以人體呼出的飛沫為載體實現人群之間的傳播，而人體飛沫在室內的傳播輸運特性受唾液滴的釋放位置，釋放速度，粒徑大小以及室內通風方式的影響[1-3]，Chen 等[4]人在分析了影響飛沫液滴擴散的眾多因素后指出室內通風方式和初始呼氣速度的影響占主導地位。當前所進行的飛沫液滴等顆粒污染物的擴散研究大多基于機械通風背景，而機械通風背景下的飛沫等顆粒污染物的擴散分析并不一定都適用于自然通風，因此有必要對自然通風背景下的飛沫液滴進行擴散分析。學校教室內人員相對密集，且大多數不具備機械通風條件，本文以自然通風教室為背景，數值分析教室進風風速對人員咳嗽飛沫液滴顆粒擴散的影響。

1 物理數學模型

1.1 物理模型

原型為東華大學第一教學樓1107 教室，桌椅為組合式桌椅，未固定可靈活擺放，具體的擺放方式為9排5 列。根據實地測量結果對教室及教室內桌椅、門窗進行尺寸取整，桌椅均簡化為實體方塊，人體模型則參考相關文獻[5]進行簡化，僅考慮人體模型的身軀和頭部。實際情形下教室開窗一側的外墻由磚墻和玻璃外墻組成，而開門一側的內墻上方位置設有常年關閉的廊窗，本文暫時不考慮廊窗和外墻上懸窗的開啟。模型及相關尺寸數據如圖1 和表1 所示。

圖1 教室模型示意圖

表1 人體及教室內相關布置物個數與尺寸

1.2 數學模型的選擇

模擬飛沫液滴擴散的CFD 方法有歐拉-歐拉法和歐拉-拉格朗日法兩種，前者將飛沫顆粒相看為連續介質，可以預測飛沫顆粒濃度分布，而后者將飛沫顆粒相看作離散相，可以對每一個粒子的運動軌跡進行跟蹤[6]。何其斌等人[7]指出由于歐拉法忽略了室內顆粒物的顆粒特性，故只適用于模擬微小的能夠較好跟隨氣流流動的中性浮力顆粒物，而拉格朗日法則主要用于預測粒子的運動軌跡和瞬態運動。本文主要研究的時咳嗽飛沫液滴的瞬態擴散過程，因此選用歐拉-拉格朗日法對咳嗽飛沫顆粒進行軌跡跟蹤研究。對于室內連續場，綜合計算速度和流場模擬精度要求，選擇計算速度適中且能很好地預測自然通風空間流場的基于RANS 方程的標準k-ε 模型進行模擬[8]。

2 邊界條件設置

本文研究的是夏季呼吸道傳染病疫情爆發期間教室前后門、4 扇平推外窗完全開啟時，自習狀態下的咳嗽飛沫液滴的瞬態擴散。設定自習狀態下教室內人數為七人，位置編號從1 到7，具體的人員分布如圖1所示。七人中僅有一人為潛在患者（潛在患者標紅），忽略其他人員咳嗽或呼吸時的液滴擴散影響，在設定的不同進風風速下分別對潛在患者處于位置1、位置3時的咳嗽液滴擴散進行分析。正常呼吸過程中產生的顆粒物液滴較少對室內影響較小，而咳嗽會產生更多的微粒或液滴，這些微?；蛞旱螘鞑ズ苓h的距離[9]，所以本文忽略潛在患者正常呼吸時產生的飛沫液滴擴散。

2.1 連續場邊界條件設置

自然通風空氣流動的動力為室內外壓力差（熱壓、風壓或二者都有），空氣流動物理過程的復雜性以及室外環境的多變性導致通風率難以預測，因此，控制自然通風以獲得所需的室內環境條件具有挑戰性，并且自然通風無法控制通過門窗口的氣流方向，室內自然通風流場的模擬難度較大[10-11]。而有關自然通風的模擬方法主要有兩種，一種是建立建筑的實際模型，根據當地的氣象資料直接給定通風口的風度大小和溫度。另一種是在建筑周圍加一個外場，擴大計算模型，在外場的入口給定當地氣象條件下對應的速度大小和溫度[12]。如果對進口風速、風向、溫度進行實時模擬，不僅計算成本巨大而且也很難實現，因此本文采用第一種自然通風模擬方法對進風邊界條件進行簡化。根據當地氣象條件假設入口進風平均速度（方向垂直于進風口），并假定進風參數在模擬過程中不發生變化，即在特定的通風模式下進行穩態模擬。

具體數值設置參照《民用建筑供暖通風與空氣調節設計規范》[13]選取，上海地區夏季室外平均風速為3.1 m/s，本文分別選取1 m/s，2 m/s 的進風風速進行研究，進風溫度和相對濕度選用夏季通風室外計算溫濕度，數值分別為31.2 和69%。進風口選用速度進口（velocity-inlet），水蒸氣質量分數計算約為0.0194，出風口（門）采用出流邊界條件（outflow），認為出口流動達到完全發展狀態。教室墻壁、天花板、地面、桌椅表面設置為絕熱壁面邊界條件，僅考慮人員顯熱忽略潛熱和輻射，人員皮膚發熱率取40 W/m2[14]。

2.2 離散相邊界條件設置

關于咳嗽過程的飛沫液滴，已有大量文獻使用不同的方法進行了研究，加之人員個體存在差異使得研究結果略有不同。孫煒等人[15]發現初始直徑大于150的飛沫液滴釋放之后很快會沉降到了地板上，Chao 等人[16]通過粒子圖像測速技術（PIV）發現健康志愿者咳嗽時平均呼氣流速為11.7 m/s，而Zhu 等人[1]所測量的咳嗽平均氣流數據為11.2 m/s。實際咳嗽過程中的呼氣速度是變化的，本文對人體以咳嗽方式產生飛沫的過程進行簡化——飛沫僅從嘴巴以水平方向咳出，飛沫速度與咳嗽呼氣速度相同，取平均流速11.7 m/s，一次咳嗽持續時間設定為0.5 s[17]，潛在患者進行間隔時間為0.5 s 的兩次咳嗽，咳嗽噴出的氣流相對濕度設定為100%，氣流與飛沫顆粒溫度均取體溫309 K，呼出氣流中的水蒸氣質量分數計算約為0.0372。

本文選取粒徑范圍為2～150 μm 的飛沫液滴為研究對象，此范圍的飛沫液滴短時間內不會迅速沉降到地面上[15]，一旦攜帶有病原體微生物則潛在危害極大。數值模擬時將飛沫液滴的物性簡化為水，密度稍微夸大設定為1003 kg/m3，蒸發溫度設定為273 K，蒸發后的液滴核粒徑設為初始液滴粒徑的32%[18]，通過計算得出蒸發性組分（水）的質量分數約為96.6%。結合Chao 等人[16]總結的咳嗽產生的不同粒徑范圍的液滴數量，依據密度估算出各個液滴組的質量，再根據飛沫液滴對應的質量百分數采用python 軟件對在2～150 μm 范圍內非均勻分布的飛沫液滴進行Rosin-Rammler 粒徑分布擬合，質量分數表和擬合圖像見表2 和圖2。

表2 一次咳嗽飛沫液滴質量分數表

圖2 咳嗽飛沫液滴粒徑分布擬合

此外，設定飛沫顆粒在門窗開口，嘴巴開口處的離散相邊界條件為“escape”；在室內人體表面、桌椅表面、地板的離散相邊界條件為“trap”，即沉降到接觸表面；而各個壁面則設為彈性壁面條件“reflect”，即顆粒到達壁面會反彈回室內空間，考慮到顆粒在壁面切向和法向上的動量損失，恢復系數均減小為0.5。

3 結果與分析

本文模擬了自然通風教室內兩種進風風速下，潛在呼吸道疾病患者分別在位置1，位置3 進行連續兩次咳嗽后釋放的特定粒徑范圍飛沫液滴的分布擴散情況。由于計算機性能和內存有限，從第一次咳嗽起共模擬了100s 的咳嗽液滴擴散情形。數值模擬首先利用SIMPLEC 算法求解穩態的自然通風室內連續相流場，之后流場中加入飛沫液滴離散相并采用PISO 算法進行非穩態的兩相相間耦合計算[19]，直至收斂。不同進風風速下的兩種情形下的咳嗽飛沫液滴分布圖如圖3 所示。

圖3 潛在患者位于1、3 位置時不同進風風速下的飛沫液滴分布圖

根據飛沫液滴分布圖可以看出，在位置1、3 處咳嗽釋放的飛沫液滴伴隨著進風風速的增大，空間殘留的飛沫液滴數量明顯減少，位置1 右前方為飛沫液滴的主要影響區域，而在位置3 釋放的飛沫液滴，由于中間兩扇外窗的間距較大，位置3 沒有直接位于氣流流動路徑上，位置3 附近存在氣流回流，飛沫液滴向位置3 左右兩側空間擴散，液滴彌漫于整個教室空間，風速一定程度上加快了飛沫向左右側的擴散。

當人員咳嗽位置確定后，進風風速是唯一變量，而進風風速對教室空間殘存液滴數量的影響，可以用特定時刻的室內空間液滴包（parcels）個數與第二次咳嗽結束后的室內空間液滴包個數的比值的變化來表征，詳見圖4。每個液滴包可以包含多個但不限于整數個飛沫液滴顆粒。從圖4 可以看出，進風風速增大可以加速1、3 位置處釋放飛沫液滴的空間數量減少速度，通過對飛沫液滴的沉降，捕獲，逃逸作用的促進來實現這一過程。與1 位置咳嗽情形不同，3 位置處釋放的飛沫液滴在100 s 后的殘存數量隨著風速增大明顯減少。由于教室空間較大，為了進一步了解咳出飛沫液滴的擴散情況，對潛在患者位于1、3 位置時不同進風風速下的坐姿呼吸平面高度（y=1.2 m）上的飛沫液滴濃度云圖進行分析，液滴濃度云圖如圖5 所示。本文對濃度云圖中的最大液滴濃度進行了縮放（由實際值調整縮小為1.00e-10kg/m3），以便于更好的顯現出飛沫影響范圍隨時間和進口風速的變化情況。

圖4 不同進風風速下的教室空間飛沫液滴數量變化情況

圖5 潛在患者位于1、3 位置時不同進風風速下坐姿呼吸高度平面上的飛沫濃度云圖

從飛沫液滴濃度云圖可以看出，當進風風速為1 m/s 時，1 位置處釋放的飛沫在呼吸平面上的濃度影響大約60 s 后很輕微，而增大風速為2 m/s 后該時間約提前到了30 s，飛沫影響區域局限于位置1 右前方。對于3 位置咳嗽情形，飛沫液滴的影響范圍更大，影響時間更長，并且飛沫對于3 位置左側附近空間呼吸高度平面區域的影響明顯強于右側，增大進風風速同樣會加快飛沫液滴濃度向3 位置左側的擴散。對于咳嗽位置周圍飛沫液滴的具體變化情況，本文在潛在患者左右方坐姿呼吸平面高度上各截取一個1200 mm×800 mm 長方形平面，分析截取平面上飛沫液滴平均濃度隨時間和進風風速的變化情況。位置3 釋放的飛沫液滴向該位置左右兩側擴散，并有可能對位置2，位置4 上的人員產生影響，在位置2，位置4 前方坐姿呼吸平面高度截取相同的長方形平面，依據是否出現平均濃度非零情形來評估飛沫影響。坐姿呼吸高度小截面上的飛沫液滴平均濃度變化折線圖及相關散點圖如圖6 所示。

圖6 坐姿呼吸高度(y=1.2 m)小截面上的飛沫液滴平均濃度變化情況

折線圖a 為潛在患者在位置1 咳嗽時兩種進風風速下的右方呼吸高度小平面上的平均液滴濃度隨模擬時間的變化情況。折線圖b，c 分別為潛在患者在位置3 咳嗽時兩種進風風速下的位置3 右方和左方呼吸高度小平面上的平均液滴濃度的變化情況，折線圖d，e 分別為1 m/s 和2 m/s 進風風速下的位置3 左方和右方呼吸高度小平面平均飛沫液滴濃度的變化情況。散點圖f 為潛在患者在3 位置咳嗽時100 s 模擬時間內飛沫液滴對位置2，位置4 上同室人員影響情況分析。

依據圖像可以看出各個呼吸高度小平面上的平均飛沫液滴濃度在1.5 s 左右到達峰值，并且峰值濃度隨進風風速增大而增大。對于在位置1 咳嗽的情形，2 m/s 進風風速下的平均液滴濃度變化較為劇烈，大約在4 s 時兩種風速下位置1 右方的平均液滴濃度開始平齊。對于在位置3 咳嗽的情形，當進風風速由1 m/s增大至2 m/s 后，初始時刻位置3 右側平均液滴濃度值總體上稍有增加，而位置3 左側平均液滴濃度則在8s 之前有大幅度提升。此外，1 m/s 進風風速下，初始時刻位置3 左右兩側平均飛沫液滴濃度差值較小，大約5s 以后左側液滴平均濃度開始一直高于右側，直至兩方濃度平齊趨向于0。2 m/s 進風風速時，初始時刻位置3 左右側平均飛沫液滴濃度差較大，左側濃度一直高于右側，大約8s 以后兩側液滴濃度開始平齊并趨向于0。位置2 和位置4 分別位于距位置3 較遠的左右兩側，根據散點圖f 可知：位置3 釋放的飛沫液滴會對位于位置2 和位置4 的同室人員產生影響，盡管影響可能很輕微。并且相較于位置2，飛沫液滴會很快擴散至位置4 呼吸高度平面，但隨著進風風速的增加，大風速下的位置2 前方的最大平均液滴濃度要高于位置4 前方的最大平均濃度，飛沫對位置2 的影響要強于位置4。

4 結論

本文在模擬過程中對進風狀況，教室熱邊界條件，人員的呼吸狀態，飛沫液滴的物性參數均進行了一定簡化。得出的主要結論如下：

1）進風風速增大可以加快1、3 位置處釋放的飛沫液滴數量的減少。但相比3 位置處釋放的飛沫液滴，模擬時間內增大進風風速對1 位置處釋放的飛沫液滴的最終殘存數量影響不大。

2）對于位置1 情形，大風速時位置1 右方呼吸高度小平面上的平均液滴濃度在4 s 之前始終比小風速下的高。對于位置3 情形，隨著進風風速增加飛沫液滴向左側的擴散效應加劇，左側平均液滴濃度在8 s之前均有大幅度提升，且同一進風風速下位置3 的左側平均濃度在初始時刻一直高于右側，大進風風速時左右兩側的平均飛沫液滴濃度差較大。

3）位置1 釋放的飛沫液滴僅對該位置右前方區域產生影響，而位置3 釋放的飛沫液滴影響范圍較廣，會對位于位置2 和位置4 的同室人員產生影響。進風風速的增加會加快位置3 釋放的飛沫液滴向位置2，位置4 呼吸平面的擴散，飛沫向位置4 呼吸平面的擴散速度相對要快些，但大進風風速下位置2 前方飛沫液滴最大平均濃度要高于位置4 前方的最大平均濃度，飛沫在大風速下對位置2 的影響較強。