含雙側分支受限空間油氣爆炸火焰行為與超壓特性大渦模擬*

劉 沖，杜 揚，梁建軍，張培理，孟 紅

（1. 陸軍勤務學院油料系，重慶 401311；2. 陸軍勤務學院勤務指揮系，重慶 401311）

汽油作為一種工業制成品，廣泛應用于交通運輸和工業清洗行業。但是，作為一種易燃易爆易揮發的危險品，在儲存、運輸及加注的全過程中，汽油時刻有著安全威脅。為有效應對安全威脅，降低人類的生命和財產損失，非常有必要針對汽油蒸氣/空氣混合氣的爆炸特征開展深入研究。此前，已有大量針對不同可燃氣體爆炸特征的研究，發現受限空間結構和形狀特征對可燃氣爆炸超壓和火焰的變化規律有較大影響。主要影響因素包括為狹長受限空間長徑比[1-2]、半受限空間開口尺寸及位置[3-5]、受限空間障礙物形狀及尺寸[6-8]和所含彎管及分支結構特征[9-11]等。特別地，彎管及分支結構是儲油洞庫坑道及輸油管路中常見的結構形式，研究表明彎管及分支結構的存在對爆炸超壓和火焰傳播速度有較強的激勵作用，并在一定程度上可與擋板式障礙物類比[9]。Zhu 等[12]發現含彎管結構受限空間可燃氣爆炸規律主要受彎管結構數量影響，而受彎管方向影響較小。Niu 等[13]對甲烷在輸氣管網中的爆炸特征進行了研究，發現爆炸超壓峰值及火焰傳播速度受連接平行管網的分支長度影響，兩分支長度不同時爆炸更為劇烈復雜。

盡管彎管和分支結構在形狀特征上有一定相似，但由于分支結構具有更復雜的通道和多個流動終點[14]，因此確有必要對含多分支結構受限空間油氣爆炸特征進行研究。Zhang 等[14]對5 種不同分支結構受限空間爆炸進行研究，發現油氣爆燃指數的增加與分支數量呈函數關系。Li 等[15]、杜揚等[16]通過實驗研究，發現分支結構可以增大汽油蒸氣/空氣混合氣爆炸超壓峰值、升壓速率及火焰傳播速度。杜揚等[17]對含雙側分支受限空間油氣爆炸進行研究，發現火焰在含有雙側分支結構的管道內會呈現“浪花狀”火焰。目前，人們已經對含分支結構受限空間油氣爆炸進行了一定的研究，但研究多集中于實驗分析和現象描述，對爆炸流場和火焰的耦合關系仍需深入研究。

本文中，在已有實驗研究[17]的基礎上，通過大渦模擬方法對含雙側分支受限空間油氣爆炸過程進行模擬，并將火焰形態、傳播速度和爆炸動態超壓等模擬結果與實驗進行驗證，分析含雙側分支受限空間油氣爆炸流場結特征、火焰傳播行為和動態超壓特性，指出“浪花狀”火焰形成的內在原因，揭示流場、火焰和超壓的耦合機理。

1 數學模型

1.1 大渦模擬控制方程

大渦模擬通過對小尺度湍流脈動進行過濾，求解計算量比直接模擬方法小，但空間分辨率比雷諾平均方法高，適合于捕捉油氣爆炸的湍流特征[18]。

大渦模擬對小尺度湍流脈動進行Favre 過濾后，利用N-S 方程對大尺度湍流脈動進行直接求解。經Favre 過濾后，大渦模擬控制方程的連續性方程、動量方程、能量方程、反應進程變量及其守恒方程為：

式中：ρ 為密度，t 為時間，ui、uj為速度分量，xi、xj為方向分量，τij為黏性應力張量分量，p 為壓力，hs為顯焓，λ 為熱導率，T 為溫度，D 為擴散系數，為歸一化化學反應速率。Yf為局部未燃燃料質量分數，Yf0為初始未燃混合氣中燃料質量分數；c 為反應進程變量，c=0 代表化學反應未進行，c=1 代表化學反應完全；符號～表示Favre 過濾量，—表示物理空間過濾量。

為使大渦模擬控制方程封閉，需建立亞網格模型。壁面自適應局部渦黏模型（WALE）能夠較好地捕捉流場從層流到湍流的轉變，且計算量較小[19]。因此，采用WALE 模型作為計算的亞網格模型。渦黏模型方程為：

式中：δij為單位張量分量；μS為亞網格黏性系數，為：

1.2 亞網格燃燒模型

湍流燃燒過程中，火焰鋒面厚度較小，為取得精確的結果，必須采用更精細的網格。但這樣會使網格量增大，運算成本陡增。因此，采用Zimont 模型[20]，對火焰鋒面進行增厚處理，確保在有限運算量下精確模擬火焰鋒面。在Zimont 模型中，湍流火焰速度公式為：

2 數值模擬

2.1 物理模型

圖1 實驗系統Fig. 1 Experimental system

以文獻[6]中的實驗數據為基礎開展數值模擬。其實驗系統如圖1 所示，由實驗管道、點火系統、壓力采集系統、高速攝影系統、油氣產生系統、油氣濃度測試系統、同步觸發系統和數據采集系統等部分組成。圖2 為實驗管道尺寸，數值模擬幾何模型與實驗管道尺寸一致。主管道及兩側分支管道內部橫截面尺寸均為100 mm×100 mm，主管道左端面密閉，中心處為點火點，右端面開口率為50%。實驗前，用聚乙烯薄膜密封右端開口，但爆炸破膜壓力很小，對爆炸發展演變影響較小，因此在數值計算過程中未考慮薄膜的影響[8]。為保證計算結果的準確性，數值計算中在管道右端開口處設置了一個700 mm×500 mm×800 mm 的外部空間模擬實驗外場，以抑制爆炸回流對管道內部的影響，如圖2 所示。沿爆炸傳播方向，在分支前、后各275 mm 的主管道處分別設置兩個壓阻型傳感器P1、P2，對應數值計算相同位置兩個壓力監測點。實驗時，管道充裝汽油體積分數1.73%的汽油-空氣混合氣，環境壓力約為1.013×105Pa，環境溫度約為300 K。高速攝影儀記錄爆炸過程火焰發展過程，幀率為500 Hz。

2.2 網格劃分

圖2 實驗管道Fig. 2 Experimental pipe

大渦模擬計算過程中，對較大尺寸渦團直接求解、對較小尺寸渦團過濾的特性決定了利用大渦模擬進行數值計算沒有“網格無關性”，越小的網格尺寸和時間步長會獲得越精確的計算結果，這與雷諾平均方法不同。但隨著網格尺寸趨近更小，計算量將趨近于直接模擬方法，因此在網格劃分上要盡可能平衡計算精度和計算成本[21]。經過調試驗證后，本文中網格劃分如圖3 所示，圖3(a)為計算域整體網格劃分，圖3(b)為分支處局部加密網格。因計算區域規則，內場和外場網格均為六面體結構網格。在內場分支結構區域加密網格，最小為2 mm；外場沿爆炸傳播方向網格逐漸增大，最大為10 mm：計算域網格總數約為490 萬個。

圖3 網格劃分Fig. 3 Computational grids

2.3 邊界條件和初始條件

為與實驗條件相對應，對含雙側分支結構受限空間充裝汽油體積分數1.73%的汽油-空氣混合氣的爆炸過程進行模擬。主管道和分支管道壁面設置為絕熱無滑移壁面且溫度為300 K，外場空間壁面設置為壓力出口且壓力為1.013×105Pa。環境溫度為300 K，流場初始速度和反應進程變量為零。在實驗管道左端面中心處設置一個半徑5 mm 的半球形區域，將其內部反應進程變量設為1，用以模擬實驗點火過程。

2.4 計算方法

計算過程中，采用二階迎風格式對對流項進行離散，采用二階中心差分格式對擴散項進行離散，采用SIMPLE 算法求解壓力和速度耦合方程組。時間步長為10-4s，每時間步迭代20 次，保證能量殘差小于10-6，動量殘差小于10-5，質量殘差小于10-6，反應進程變量殘差小于10-4。所有計算在Think Station P700 工作站完成，耗時約100 h。

3 結果分析

油氣爆炸關鍵特征通常包括火焰傳播行為和超壓變化規律，為驗證數值模擬結果的準確性，選取火焰傳播形態、火焰鋒面傳播速度和超壓曲線三方面進行數值模擬，并和文獻[17]實驗數據對比。

3.1 火焰傳播行為

圖4 為火焰形態的實驗和數值模擬結果。在初始階段（32 ms 前），火焰受到管道側壁的約束沿軸向拉伸，以典型的“指形”火焰傳播，且此時火焰鋒面未受到分支結構的擾動作用，火焰鋒面包絡面清晰光滑，為層流狀。當火焰傳播至分支管道附近時，因壁面約束消失和分支管道尖角的障礙物作用，致使流場波動劇烈，火焰從層流狀向湍流狀轉捩，火焰鋒面褶皺彎曲。火焰傳播至分支管道（54～64 ms）時，火焰鋒面形成了三路傳播路徑，即火焰繼續沿主管道向末端傳播，同時在兩側分支管道中對稱傳播。在分支管道中，火焰鋒面因受到右側壁面的約束，在傳播中逐漸向左側壁面彎曲，并形成了實驗中對稱的“浪花狀”火焰。并且，隨著反應的不斷進行，狹長火焰鋒面兩側大量未燃氣被同時引燃，火焰鋒面嚴重變形并形成明顯褶皺。

圖4 實驗和數值模擬的火焰形態Fig. 4 Experimental and simulated flame structures

由圖4 可以看出，火焰在主管道中傳播時，數值模擬能夠直觀準確地反映實驗中的火焰形態變化。當火焰傳播進入分支管道時，數值模擬能夠較準確地反映實驗中“浪花狀”火焰的輪廓和鋒面位置，但在細節上和實驗仍存在一定差別。究其原因有兩點：一方面，高速攝影儀在拍攝較大尺寸物體時像素數會呈指數降低，實驗照片對火焰的細節捕捉不夠；另一方面，數值模擬采用了Zimont 模型計算火焰傳播鋒面，而該模型為了更好地捕捉到火焰鋒面位置，對火焰鋒面進行了增厚處理，導致在60 和62 ms 時對分支管道中火焰鋒面末端彎曲細節表達不精準。鑒于本文的研究目的，是為了分析火焰與流場耦合規律，因此可以忽略。

在火焰形態對比基礎上，對火焰位置和傳播速度進行了對比，如圖5～6 所示。由圖5 可以發現，數值模擬的火焰位置能夠較好地反應火焰傳播位置隨時間的變化，最大誤差為14.436%，發生在42 ms處。這是因為，此時火焰鋒面傳播至分支管道附近，爆炸流場復雜，火焰鋒面結構出現嚴重潰散，擾動隨機性較大，導致誤差增大。由圖6 可以發現，數值模擬的火焰速度與實驗的變化規律整體一致，在爆炸初期（28 ms 前），因參與反應的油氣較少，反應速率較小，火焰傳播速度緩慢升高；隨著化學反應持續進行（實驗28～40 ms，數值模擬28～45 ms），火焰鋒面向末端傳播過程中沿軸向不斷拉伸，火焰表面積增大，參與反應的單位時間油氣增多，化學反應速率增大，進而火焰傳播速度迅速增大；當火焰鋒面傳播到分支管道附近時，初期分支管道的泄壓作用較大，火焰傳播速度增長趨勢出現停滯，其后隨著分支管道內油氣被引燃，火焰鋒面表面積呈指數增加，火焰傳播速度繼續增大。數值模擬誤差最大為24%，發生在45 ms 處。

圖5 實驗和數值模擬的火焰位置Fig. 5 Experimental and simulated flame front locations

圖6 實驗和數值模擬的火焰速度Fig. 6 Experimental and simulated flame speeds

3.2 動態超壓變化規律

在壓力傳感器P1 處實驗和數值模擬的動態超壓曲線，如圖7 所示。表1 為實驗和數值模擬的3 個典型超壓峰值的對比。可見，兩者超壓曲線的規律基本一致，分別在65 和68 ms 達到爆炸超壓峰值，峰值誤差為8.6%。但在數值模擬中，未考慮薄膜破裂時的破膜壓力對管道內壓力場的影響，因此在爆炸超壓上升階段，實驗存在3 個峰值，而數值模擬只存在2 個峰值。

表1 實驗和數值模擬的典型超壓峰值Table 1 Experimental and simulated typical overpressure peaks

在不考慮破膜壓力的情況下，爆炸超壓升壓過程可以劃分為4 個階段，以數值模擬得到的曲線為例：第Ⅰ階段為緩慢升壓階段（0～20 ms），該階段從模擬點火開始至火焰拉伸呈“指形”終止，在該階段爆炸集中于點火源附近，參與反應未燃氣量少，反應釋能小，超壓隨之上升緩慢；第Ⅱ階段為超壓上升階段（20～50 ms），隨著爆炸進行，參與反應的未燃氣逐漸增多，單位時間釋能迅速增大，超壓開始明顯上升；第Ⅲ階段為超壓平衡階段（50～56 ms），主管道泄流效應和分支管道二次泄壓作用弱化了內場超壓上升，超壓上升速率減慢；第Ⅳ階段為快速升壓階段（56～65 ms），分支管道增大了燃燒面積，反應速率迅速上升，超壓再一次快速上升，直至可燃氣消耗殆盡轉入下降階段。

在小尺度狹長空間發生油氣泄爆時，通常管道內壓力波動以靜壓為主，動壓對超壓的影響非常小，因此空間內部各點的超壓基本一致。而對于含雙側分支的受限空間，當壓力波傳播經過分支附近時，內部空間橫截面的變化必然影響壓力波的傳播，因此在分支管道前后測點的超壓一定會有所差異。為此，對比了分支管道前后測點P1 和P2 處的超壓曲線，如圖8 所示。由圖8 可以看出，P1 和P2 處的超壓曲線演變規律基本一致，僅在局部有數值大小等方面的差異，具體有3 點：（1）在緩慢升壓階段和超壓上升階段，P2 處超壓上升滯后于P1；（2）P1 處最大超壓峰值大于P2，形成時刻略滯后于P2；（3）P2 處的超壓振蕩比P1 劇烈。P1 和P2 處的超壓特性差異說明，分支管道導致主管道橫截面突擴，在分支管道附近誘導稀疏波導致壓力波衰減，分支管道對主管道起到了泄壓作用，并且P2 位置更靠近泄爆口，受泄流效應的影響也更大。

圖7 P1 處實驗和數值模擬的爆炸超壓曲線Fig. 7 Experimental and simulated overpressure curves at point P1

3.3 流場-火焰-超壓耦合

圖8 P1 和P2 處數值模擬的爆炸超壓曲線Fig. 8 Simulated overpressure curves at points P1 and P2

由文獻[22]，在爆炸發展過程中，流場變化會影響火焰傳播行為，火焰傳播行為的改變亦會對爆炸流場產生影響，火焰傳播和湍流流場間存在一定的正反饋耦合機制。此外，火焰結構在一定程度上表征反應速率的快慢，進而影響爆炸超壓的波動。

圖9 為火焰傳播形態和流場速度矢量圖，白色箭頭代表速度矢量，顏色越白速度數值越大，速度矢量所在平面為管道中軸面（y=0.05 m）。由于重點關注分支管道附近的火焰和流場形態，因此火焰傳播至分支管道前（20、32 ms）的圖片包含點火端至分支管道；而為了清晰顯示火焰在分支管道中傳播形成“浪花狀”的過程，這個階段（48、60、64 ms）的圖片僅包含分支管道區域。

在20 ms 時，火焰以層流狀態傳播，火焰鋒面前端附近速度矢量沿半球面指向未燃區，主管道未燃區氣體在火焰膨脹作用下沿水平方向向管道末端流動。當流動經過分支管道時，兩側壁面對未燃氣的約束消失，未燃氣對稱地向兩側分支管道膨脹，在主管道與分支管道交界處建立了層流流場，此時分支管道內部流場速度趨近于零。隨著反應進行，在32 ms 時火焰以“指形”傳播，火焰面積增大，主管道軸向速度場迅速增強，由于受到分支管道尖角的擾動流場在分支管道口附近向湍流轉捩，且因分支管道右側壁面的約束，在主管道與分支管道交界處形成了對稱的渦旋結構，上側旋轉方向為逆時針，下側為順時針，此外分支管道內部流場速度逐漸增大。當反應進行到48 ms、火焰傳播至分支管道口時，火焰鋒面在復雜湍流流場作用下顯著變形，而早期建立的對稱渦旋結構在火焰沖擊作用下向分支管道內部發展且旋轉方向維持不變，但渦旋結構沿分支管道軸向拉伸成為較大的橢圓形，此時火焰面積較小且受到較強的泄壓作用，導致火焰傳播速度上升趨勢減緩。此后，火焰傳播進入分支管道，在漩渦結構作用下火焰鋒面沿分支管道右側壁面拉伸，而后向左側壁面彎曲形成60 ms 時的“浪花狀”火焰。與此同時，火焰鋒面前側未燃氣不斷被引燃，火焰膨脹促進分支管道內部速度矢量指向主管道，分支管道內未燃氣在此流場作用下進入主管道，導致主管道與分支管道交界處形成新的燃燒面，如64 ms 時。

為了便于分析火焰形態和動態超壓的耦合關系，將P1 處數值模擬的動態超壓曲線和典型時刻火焰形態結構一起呈現，如圖10 所示。根據前述分析，爆炸超壓變化過程可劃分為4 個升壓階段和1 個振蕩下降階段。第Ⅰ階段，在點火初期，只有點火頭周圍少量油氣參與化學反應，而后火焰以半球形傳播，參與反應的可燃氣有限，爆炸緩慢發展，超壓上升趨勢不明顯。第Ⅱ階段，爆炸進入快速升壓階段，初期火焰以“指形”結構傳播（32 ms），火焰面積隨著鋒面拉伸不斷增大，參與反應的油氣增多，超壓逐漸上升，當火焰傳播至分支管道附近時，因流場復雜導致的火焰鋒面褶皺嚴重進一步導致反應更加劇烈，超壓上升趨勢更明顯。第Ⅲ階段，火焰傳播進入分支管道形成“偽受限空間”燃燒（54、56 ms），初始火焰傳播速度較慢，火焰面積變化不明顯，參與反應的油氣量沒有顯著增加，而火焰鋒面靠近末端泄爆口，主管道泄壓作用增強，導致這個階段升壓速率下降。第Ⅳ階段，當火焰在分支管道中渦旋結構的作用下充分發展時，火焰鋒面不規則潰散嚴重，同時參與反應的油氣量也大幅增多，爆炸劇烈進行，升壓速率較前3 階段都大。振蕩下降階段，火焰傳播到64 ms 時，由于主管道末端的泄壓效應，超壓迅速下降，而后在壓力波來回傳播和管道內殘余油氣持續爆炸的兩方面作用下形成超壓振蕩。

圖10 數值模擬的超壓曲線和火焰傳播形態Fig. 10 Simulated overpressure curve and flame propagation structure

4 結 論

利用大渦模擬WALE 模型，模擬了雙側分支結構受限空間油氣泄壓爆炸的超壓特性和火焰行為，與實驗數據對比了火焰傳播形態、火焰傳播速度和動態超壓曲線，驗證了所選取模型對含多分支結構受限空間泄壓爆炸計算的有效性。然后，運用該模型對雙側分支結構受限空間油氣泄壓爆炸的空間流場結構、火焰行為和超壓特性進行了研究。得出以下結論：

（1）火焰傳播進入分支管道前，火焰行為會影響分支管道內部流場，在主管道和分支管道交界處會產生方向相反的對稱渦旋結構，隨著火焰進一步傳播，渦旋結構會向分支結構內部發展，建立旋轉流場；

（2）當火焰傳播進入分支管道后，分支管道內部前期已建立流場決定了火焰的形態，火焰鋒面在渦旋結構作用下形成“浪花狀”火焰，此后火焰和流場相互影響，流場向湍流轉捩，火焰鋒面潰散；

（3）爆炸超壓升壓過程可劃分為4 個階段，受到火焰鋒面面積、主管道泄流效應和分支管道二次泄壓共同作用，爆炸流場、火焰行為和動態超壓具有顯著耦合性。