發煙罐煙霧質量濃度時空分布特性

楊尚賢，陳慧敏，高麗娟，馬超，齊斌，鄧甲昊

(北京理工大學 機電動態控制重點實驗室,北京 100081)

0 引言

煙霧能夠在一定時間和空間范圍內對可見光、激光、毫米波等形成有效干擾和衰減。同時，人工煙霧釋放操作簡單、成本低廉，在現代戰爭領域得到了迅猛發展與廣泛應用。研究煙霧的釋放、擴散規律和濃度分布狀態，為獲取激光在煙霧環境下的傳輸特性，提高對地激光類武器系統的抗煙霧干擾性能，具有重要意義[1-3]。

發煙罐的作用原理是罐中發煙劑(以磷基發煙劑為代表)與氧氣作用而燃燒，生成具有高分散度的氧化物粒子并擴散形成煙霧。煙霧較長時間懸浮在大氣中，實現遮蔽、迷盲或干擾制導的目的[4]。為了得到煙霧的擴散規律、探索煙霧質量濃度的時空分布特性，近年來學者們做了大量研究。文獻[5]采用數值模擬和野外測試相結合的方法，將煙幕初始云團的爆炸分散過程分為內部壓力衰減為0之前和之后兩個階段，得到兩個階段煙幕初始云團的最大半徑。文獻[6]應用計算流體力學仿真軟件Fluent研究了風速和風向對爆炸式催淚彈煙霧擴散的影響。文獻[7]采用瞬時點源模擬爆炸煙霧，基于高斯煙團模式進行建模仿真，得到了煙霧濃度分布規律。文獻[8]基于解析幾何分析構建了垂直煙幕空間尺度計算模型，并基于高斯煙羽模式，應用牛頓下山法構建了單點釋放垂直煙幕計算模型。文獻[9]分析了發煙彈擴爆和煙霧在大氣中擴散兩個過程，提出了基于發煙彈擴爆后相關參數及氣象條件來估算煙霧濃度的方法。文獻[10]以湍流大渦算法對設置擋煙垂壁的防護工程煙霧擴散情況進行了仿真研究。文獻[11]對地下娛樂建筑發生火災時的煙霧流動進行模擬，分析了煙霧流動狀態并與機械排風進行了比較。文獻[12]基于氣體- 固體兩相流噴射的機理，數值模擬了微重力、高真空環境下的空間煙霧擴散。

以上文獻應用不同方法針對煙霧擴散過程進行了研究，但仍存在以下問題：1)沒有對煙霧的擴散進行分階段、分時刻的詳細研究，對煙霧擴散時煙團形狀的變化過程沒有展開描述；2)基于經驗公式構建的煙霧擴散模型，以及基于煙霧釋放參數和氣象條件估計的煙霧空間質量濃度分布，數據誤差較大；3)建筑通道對煙霧的擴散有阻礙作用，與發煙罐煙霧在自由空間中的擴散有明顯差異。另外，大部分文獻沒有分析不同環境因素對煙霧擴散的影響規律。

為得到不同因素下煙霧擴散的時空分布特性，本文基于計算流體力學，應用Fluent軟件對一種小型發煙罐煙霧的擴散過程進行仿真模擬。該方法能夠直觀地呈現計算結果，從而可以對煙霧的流動結構進行細致的分析。通過分析擴散時間、風速、質量流率等因素對煙霧形狀、煙霧質量濃度分布狀態的影響規律，為人工煙霧釋放提供理論指導，并為進一步研究激光在煙霧中的傳輸特性奠定基礎。

1 煙霧擴散數值模型

煙霧擴散過程屬于典型的離散相擴散模型，建模時要同時考慮氣流流動和煙霧顆粒的擴散運動。

1.1 氣流流動數值模型

煙霧在上升過程中受到迎面空氣的作用(包括空氣介質的阻力、能量的交換)與上升流本身的擴散作用。設定流場速度u在3個坐標軸方向的分量為ux、uy、uz，溫度場參數為T，流場壓力為p.煙霧擴散過程遵循質量守恒定律、能量守恒定律和動量守恒定律。另外，煙霧的流動狀態屬于湍流流動，工程中常用雷諾時均方程來描述湍流狀態，基本思想是通過k-ε(k為湍流動能；ε為耗散率)雙方程模型表示時均方程中流體的瞬時脈動量[13]。相關的數學模型有如下5個。

質量守恒方程：

(1)

式中：ρ為空氣密度(kg/m3)；xi為湍流模型的張量表示形式，i=x，y，z；ui(i=x，y，z)為流場速度在x軸、y軸、z軸方向上的分量；源項Sm為持續流入的氣體質量(kg)。

能量守恒方程：

(2)

式中：cp為空氣的定壓比熱容(J/(kg·℃))；λ為氣體的導熱系數(W/(m·℃))；源項ST為持續流入的氣體熱量(J)。

動量守恒方程：

(3)

式中：uj為流場速度在x軸、y軸、z軸3個方向的分量，j=x，y，z；xj為湍流模型的張量表示形式；μ為層流黏度系數(Pa·s)；μt為湍流黏度系數(Pa·s)，μt=ρCμk2/ε，Cμ=0.09.

k方程：

(4)

ε方程：

(5)

式中：C1ε=1.44；C2ε=1.9；σε=1.3.

1.2 煙霧顆粒數值模型

應用Fluent軟件中嵌套的離散相模型模擬煙霧顆粒的擴散運動。離散相模型是遵循歐拉-拉格朗日法的兩相流數值模型，該模型將流體視為連續相，在拉格朗日坐標系中對顆粒作用力微分方程進行積分來求解顆粒的運動軌道。在笛卡爾坐標系中，煙霧顆粒的作用力平衡方程和軌跡方程(x軸方向)描述如下：

(6)

(7)

式中：u1、up分別為x軸方向流體相速度和顆粒速度(m/s)；FD(ul-up)為顆粒的單位質量曳力函數；ρp為顆粒堆積密度(kg/m3)；gx(ρp-ρ)/ρp為顆粒的單位質量重力與浮力之差；Fx為單位質量的其他作用力。對(6)式在離散的時間步長上逐步積分，即得到顆粒軌道上每一個位置上的顆粒速度；對(7)式沿每個坐標軸方向積分求解，即得到顆粒相的軌跡。

2 物理模型和基本參數設置

數值求解前建立物理模型并對數值模擬的參數進行設置，物理模型和參數設置的準確性決定了數值模擬結果的可靠性。

2.1 物理模型與網格劃分

煙霧生成過程是發煙材料燃燒的過程，而燃燒是復雜的化學反應。本文重點研究發煙罐煙霧在自由空間的擴散與質量濃度分布狀態，對發煙材料的燃燒過程不做分析。因此，仿真模型簡化為帶有一定初始速度的高溫煙霧從發煙罐底部涌出，并在自由空間擴散。構建的發煙罐物理模型和計算域如圖1所示。

圖1 發煙罐物理模型和計算域Fig.1 Physical model and computational domain of smoke pot

由圖1可見，發煙罐模型為直徑70 mm、高度400 mm的圓柱，上端開口，下端封閉并置于地面上。計算域為10 m×10 m×10 m的立方體。參考坐標系Oxyz如圖1所示，發煙罐底部中心為坐標原點O.

利用ICEM軟件對上述幾何模型進行四面體網格劃分，兼顧仿真精度和收斂速度。由于發煙罐表面附近煙霧的速度、溫度梯度大，對其網格進行了加密；遠離發煙罐的區域對計算精度要求不高，網格以一定增長率逐漸變大，以減小計算量。另外，在Fluent軟件中將四面體網格轉化為多面體網格，并對網格進行光滑操作[14]。計算域網格劃分結果如圖2所示，最終計算域網格數量為28 627個，其中網格質量優于0.71和0.66的網格分別占95%和99%.

圖2 計算域網格劃分結果Fig.2 Partitioning result of computational domain mesh

2.2 基本參數設置

發煙材料燃燒后生成煙霧顆粒，其粒徑分布是煙霧的基本特性參數之一。將煙霧顆粒尺寸劃分成一系列尺寸帶，每個尺寸帶中煙霧顆粒占顆粒總量的百分數就是煙霧的粒徑分布。采集小型發煙罐燃燒釋放后的煙霧顆粒樣品，使用日本HORIBA公司生產的LA-950型激光散射粒度分布分析儀對粒徑進行測試分析，得到煙霧顆粒的最小粒徑為3 μm，最大粒徑為26 μm，平均粒徑為9 μm.煙霧的粒徑分布如圖3所示。

圖3 煙霧顆粒的粒徑分布Fig.3 Particle size distribution of smoke particles

對煙霧粒徑進行最小二乘回歸分析和顯著性檢驗，得到煙霧顆粒的粒徑服從羅辛- 拉姆勒(Rosin-Rammler)分布。羅辛- 拉姆勒分布又稱R-R分布，是最常用的粉塵顆粒分布函數[15-16]，其表達式為

F(d)=1-exp[-βdn],

(8)

式中：β為特征參數；d為粒徑尺寸；n為分布指數；F(d)為粒徑小于d的累積百分比。根據煙霧粒徑分布數據，可求解出n=3.45，β=3.09×10-4.

發煙材料的成分主要為磷與硝酸鉀、硝酸鋇、二氧化錳、三氧化二鐵、高錳酸鉀等構成的混合物，其燃燒釋放的熱量使得煙霧的初始溫度可達到400～500 ℃[17].設定仿真條件如表1所示。

表1 仿真條件設置Tab.1 Simulation conditions

3 仿真結果和分析

3.1 無風情況下煙霧的擴散分析

首先分析無風情況下煙霧釋放后的擴散狀態，選擇瞬態計算模型進行仿真。將發煙罐煙霧的擴散過程分成兩個階段：0～10 s，煙霧持續從發煙罐底部噴涌；10～100 s，停止噴涌，煙霧在自由空間中彌散。

沿z=0 m截面觀察不同時刻的煙霧狀態，煙霧質量濃度在0～10 s的分布如圖4所示。由圖4可以看出，隨著煙霧的持續噴涌，煙羽體積逐漸增大，外形呈傘狀。這是因為除迎面空氣阻力外，煙霧沒有受到障礙物和風力等干擾，煙羽中軸周圍壓力保持平衡，煙羽中軸煙霧質量濃度明顯比周圍質量濃度高。

圖4 煙霧質量濃度在0～10 s的分布Fig.4 Distributions of smoke mass concentration at 0-10 s

為了研究煙霧溫度對煙霧擴散速度的影響規律，分析3 s、5 s、8 s、10 s 4個典型時刻的空間溫度場分布和垂直方向(y軸方向)空間速度場分布，結果分別如圖5、圖6所示。圖5、圖6中，垂直方向空間速度場表征了煙霧擴散速度，煙霧上升速度越快，擴散速度越快。

由圖5、圖6可以看出，高溫煙霧從發煙罐出口噴出后與冷空氣進行熱量的傳遞，溫度降低。對比同時刻的煙霧質量濃度分布和溫度場分布可以發現，高溫煙霧主要集中于煙霧中軸，而從中軸附近擴散開后溫度基本下降到與環境溫度相同。發煙罐出口正上方，持續的煙霧供給為煙羽中軸部位提供熱量，因此中軸部位的煙羽溫度高于環境溫度并能夠以較快的速度向上噴涌。煙霧從中軸附近擴散開后不能得到持續的高溫煙霧供給，因此降溫較快，上升速度也明顯減慢。

圖5 典型時刻的空間溫度場分布Fig.5 Spatial temperature field distribution at the characteristic moment

圖6 典型時刻的空間速度場分布(垂直方向)Fig.6 Distribution spatial velocity field in vertical direction at characteristic typical moment

取樣分析3 s、5 s、8 s、10 s 4個典型時刻煙羽中軸溫度、上升速度分別與高度的變化關系，繪制曲線分別如圖7、圖8所示。由圖7可以看出：4個時刻的煙羽中軸溫度隨高度增加先迅速下降，然后下降速度趨緩，最后達到穩定狀態；4個時刻的煙羽長度依次遞增，造成4個時刻煙羽中軸溫度隨高度的下降趨勢依次放緩。對比圖7、圖8可以看出，煙羽中軸溫度、上升速度隨高度的變化趨勢基本同步。

圖7 典型時刻煙羽中軸溫度與高度的變化關系Fig.7 Relation between temperature and height in the middle axis of plume at a typical moment

圖8 典型時刻煙羽中軸上升速度與高度的變化關系Fig.8 Relation between rising speed and height in the middle axis of plume at typical moment

對每個取樣點的溫度和對應的上升速度進行曲線擬合，得到相關性曲線如圖9所示。由圖9可以看出，煙霧溫度與上升速度有明顯的線性關系。

圖9 溫度與上升速度的相關性曲線Fig.9 Correlation curve between temperature and rising speed

煙霧質量濃度在10～100 s的分布如圖10所示。由圖10可以看出，煙霧停止噴涌后，煙羽持續擴散并呈團狀。煙團中軸部分上升速度仍然較快，但是沒有了煙霧補充，導致中間煙霧質量濃度低而周圍煙霧質量濃度高，周圍煙霧圍繞中軸形成一個明顯的煙環。隨著時間的延續，煙團體積持續增大并呈云狀。整體來看，煙云上升速度呈逐漸放緩的趨勢。

圖10 煙霧質量濃度在10～100 s的分布Fig.10 Distributions of smoke mass concentration at 10-100 s

3.2 風速對煙霧擴散的影響

發煙罐在戶外開放空間使用，環境中的風極易對煙霧的擴散狀態造成影響，進而影響發煙罐的遮蔽效果。對不同風速下的煙霧分布狀態進行仿真分析，分別設定風速為0.5 m/s、1.0 m/s、1.5 m/s、2.0 m/s、2.5 m/s、3.0 m/s、4.0 m/s、5.0 m/s、10.0 m/s、15.0 m/s，風流方向沿x軸正方向。煙霧噴涌時間為10 s，仿真總時間為10 s.取不同風速下第5 s和第10 s典型時刻的煙霧質量濃度分布狀態進行對比，沿z=0 m截面觀察，結果如圖11、圖12所示。

圖11 不同風速下5 s時刻的煙霧質量濃度分布Fig.11 Distributions of smoke mass concentration at 5 s at different wind speeds

圖12 不同風速下10 s時刻的煙霧質量濃度分布Fig.12 Distributions of smoke mass concentration at 10 s at different wind speeds

由圖11和圖12可以看到，環境中的風使煙羽發生了傾斜，且隨著風速的增加煙羽傾斜幅度增加。有利于煙霧較大范圍彌漫于近地面附近，有效提高對地面保護目標的遮蔽效果。當風速小于2.5 m/s時，10 s時刻的煙羽明顯長于5 s時刻的煙羽，且擴散范圍更大，而當風速超過2.5 m/s后則差異不大，即5 s時刻煙霧已經達到了穩定擴散狀態。這是因為風加快了煙霧的擴散，使煙霧在較短時間內達到穩定擴散狀態，快速增大遮蔽范圍。但是，當風速增加到一定程度后又會造成煙霧被急劇吹散，煙霧質量濃度太低，遮蔽效果下降。對比風速為4 m/s、5 m/s、10 m/s、15 m/s時的煙霧質量濃度分布，可以發現煙霧整體上的空間質量濃度依次降低。風速在2.0 m/s、2.5 m/s、3.0 m/s時，發煙罐釋放的煙霧能在較短時間內達到穩定擴散狀態，且在穩定擴散時具有較大的遮蔽范圍和較高的空間濃度。因此，初步認為發煙罐的最佳使用風速在2～3 m/s附近。

在計算域空間中選擇典型位置進行特征點采樣，分析煙霧質量濃度與風速之間的關系。在z=0 m截面上沿煙霧擴散方向選取P(1.0 m，1.0 m，0 m)、P(2.0 m，1.0 m，0 m)、P(3.0 m，1.0 m，0 m)、P(4.0 m，1.0 m，0 m)4個采樣點，得到5 s時刻和10 s時刻質量濃度與風速的關系曲線，如圖13所示。

圖13 不同時刻4個采樣點的質量濃度與風速關系Fig.13 Relationship between mass concentration and wind speed at 4 sampling points at different times

總體來看，隨著風速增大各采樣點的濃度值先增大后減小。5 s時每個位置的質量濃度峰值，主要集中在風速為1～2 m/s時。10 s時每個位置的質量濃度峰值，主要集中在風速為2～3 m/s時。這也從數值上驗證了發煙罐的最佳使用風速在2～3 m/s附近的結論，此風速下發煙罐煙霧能在較短的時間內達到較大的遮蔽范圍和較高的空間質量濃度。

另外，對比風速為3 m/s、4 m/s、5 m/s、10 m/s、15 m/s 5種情況下煙霧已經達到穩定擴散狀態的采樣點質量濃度與風速的關系，為負相關，近似呈反比關系。這一結果符合文獻[8]引用的高斯煙羽模式中下風向任一點的煙霧濃度與風速大小呈反比的表述。

3.3 質量流率對煙霧擴散的影響

z=0 m截面不同質量流率下10 s時刻的煙霧質量濃度分布狀態如圖14所示。由圖14可以看到，隨著質量流率的增加，煙霧的空間質量濃度明顯增加。選取z=0 m和z=1.0 m兩個截面，在每個截面上選取6個采樣點，分析發煙罐質量流率變化對煙霧空間質量濃度造成的影響。

圖14 不同質量流率下煙霧質量濃度分布狀態Fig.14 Smoke mass concentration distributions at different mass flow rates

為與發煙罐的使用場景相結合，采樣點選在計算域內煙霧擴散范圍較大的區域。如圖15所示分別是z=0 m和z=1.0 m截面上(4.0 m，1.0 m，0 m)、(4.0 m，2.0 m，0 m)、(4.5 m，2.0 m，0 m)、(4.5 m，1.0 m，0 m)、(4.0 m，1.5 m，1.0 m)、(4.5 m，1.5 m，0 m)6個采樣點質量濃度值隨質量流率的變化曲線。由圖15可以明顯看出，采樣點的煙霧質量濃度值隨發煙罐質量流率的增加而線性增加。由此可以認為在1～20 mg/s的質量流率范圍內，發煙罐煙霧的空間質量濃度與質量流率呈線性關系。

圖15 不同采樣點的質量濃度隨質量流率的變化曲線Fig.15 Change curves of mass concentration at sampling point with mass flow rate

4 結論

基于計算流體力學的離散相擴散理論，本文提出一種發煙罐煙霧擴散仿真方法。應用Fluent仿真軟件分析了煙霧在噴涌和彌散階段的外形形態，溫度與上升速度的關系以及風速、質量流率等參數對煙霧擴散速度、空間質量濃度分布的影響規律。得出主要結論如下：

1) 無風時，煙霧在噴涌階段呈傘狀，停止噴涌煙霧呈團狀并存在一個質量濃度較高的煙環。煙霧從發煙罐出口噴出后，隨著高度的增加，溫度先迅速下降，然后下降速度趨緩，最后達到穩定狀態。煙霧溫度與上升速度有明顯的線性關系。

2) 一定速度的風使煙羽發生傾斜，有利于煙霧較大范圍彌漫于近地面附近，提高對地面保護目標的遮蔽效果。同時，風加快了煙霧擴散速度，仿真得到發煙罐最佳使用風速在2～3 m/s附近。此風速下發煙罐既能快速形成空間遮蔽，又具有較大的遮蔽范圍。

3) 對不同采樣點的煙霧質量濃度進行分析，發現在1～20 mg/s的質量流率范圍內，煙霧的空間質量濃度與質量流率呈線性關系。

本文仿真過程沒有考慮空氣濕度，忽略了煙霧顆粒在大氣中的吸濕性。另外，多發煙罐協作部署下的煙霧擴散規律和環境因素對遮蔽效能的影響是將來的研究方向。