泄壓條件下氣體爆炸超壓計算模型研究

袁丹燦，馬 斌，杜坤澤

(武警學院 研究生隊，河北 廊坊 065000)

泄壓條件下氣體爆炸超壓計算模型研究

袁丹燦，馬 斌，杜坤澤

(武警學院 研究生隊，河北 廊坊 065000)

基于封閉空間壓力增長理論和泄放質量比例的壓力衰減假設，建立泄壓條件下氣體爆炸超壓變化的微分方程模型，并通過固定時間差分離散將其轉化為差分方程模型，給出差分方程的求解步驟。對模型求解進行時間差分的敏感性分析，得到時間差分的大小應控制在0.001 s以下。根據Bartknecht在1 m3且有0.36 m2泄壓口的正方體空間內進行的甲烷爆炸系列試驗進行求解對比，分析了模型的適用條件、誤差因素和湍流因子的誤差機理。分析表明：對于體積為1 m3且泄壓面積為0.36 m2空間內的甲烷爆炸，在泄壓構件相對靜態觸發壓力介于20～200 kPa條件下模型適用；應用于其他條件時，應當考慮通過燃爆指數折算的湍流因子的適用范圍。

氣體爆炸；防爆泄壓；超壓模型；湍流因子

0 引言

泄壓是降低爆炸破壞最有效的手段，通過布置泄壓構件可以在爆炸全面發展之前及時將壓力泄放，從而降低爆炸可能產生的最大超壓，進而防止對結構造成破壞。工程上常根據結構強度確定泄壓條件下空間內部允許最大超壓，并以此為根據指導泄壓設計[1-3]。針對可燃氣體，Swift[4]提出泄壓面積應與允許最大超壓的平方根成反比，Bartknecht[5]則提出了該壓力與氣體燃爆指數共同決定泄壓面積的非線性經驗模型。部分學者通過研究壓力變化規律，分析壓力產生機理以確定泄壓條件下的最大超壓：Zalosh等[6-8]認為無量綱的壓力可用無量綱泄壓面積的函數關系表征，Hernandez[9]等基于已燃和未燃氣體組分變化建立了微分方程計算泄壓條件下封閉空間內壓力的演化過程，Chyzy[10]等用正弦關系描述壓力變化過程進而簡化了求解超壓的微分方程模型。本文基于一定的簡化假設，建立泄壓條件下氣體爆炸超壓的微分方程模型及求解的差分方程模型，并對模型的適用性與局限性進行探討。

1 泄壓條件下氣體爆炸超壓微分方程模型

1.1 泄壓口未開啟下的壓力增長規律

將空間內氣體視為整體，假設壓力分布均勻，在泄壓構件開啟之前空間封閉。空間內部壓力P的增長率與燃燒半徑rb有關。在長徑比L/D<2的條件下，可以認為火焰面以球面的形式從點火源向各個方向等效傳播，有

式中，t是時間，s；Su是氣體的基本燃燒速度，m·s-1；χ是湍流因子；a是空間的當量半徑，m；Pm是封閉空間內氣體爆炸產生的最大壓力，kPa；P0是初始的環境壓力，kPa；γ是氣體絕熱指數。

燃燒半徑的增長又與爆炸壓力有關

根據式(1)和式(2)，給定初始條件ε=(rb/a)→0即可求解空間壓力的變化過程。

式(1)中χ的變化規律相對復雜，與物質類型和壓力的增長過程等諸多因素相關。對于氣體的類型，可以通過立方體定律表述

式中，KG是氣體燃爆指數，bar·m·s-1。在絕熱假設下，壓力的增長不會出現拐點，因此(dP/dt)max應在P=Pmax處達到，有

故可以據此倒推估算湍流因子

將式(5)所得湍流因子視為常量，代入式(1)、式(2)求解壓力P，即可得到空間內壓力變化規律。

1.2 泄壓口開啟下的壓力變化規律

忽略泄壓構件的動作規律，假設達到其靜態觸發壓力Pstat時泄壓口完全打開，其面積為Av，m2。此時空間內氣體由于超壓將通過泄壓口向外泄放。判斷氣流的流態：若(P/P0)<[(γ+1)/2]γ /(γ-1)，為亞聲速流動，其泄放速率為

若(P/P0)≥[(γ+1)/2]γ /(γ -1)，氣流則為超聲速流動，其泄放速率為

式(6)和式(7)中，mv是泄放的氣體質量，kg；Cd是泄放系數，由泄壓構件的形式決定；ρ是空間內氣相平均密度，kg·m-3，隨壓力變化關系為

由于氣體的泄放，空間內氣相物質的質量減少，進而引起壓力增長速率衰減。假設泄壓作用下，空間內實際壓力Pv可以用泄放的質量比例系數fvent來描述，則

此時空間內的壓力增長速率也應由Pv決定，則式(1)、式(2)的P應相應作Pv，并以此求解壓力變化規律，得到的峰值壓力Pv,max對應的相對壓力Pv,max-P0就是泄壓條件下的最大超壓。

2 固定時間離散的差分方程模型

2.1 差分離散模型的建立

為便于求解，將微分方程模型基于固定時間差分Δt進行離散。此時P與Pv可統一，差分關系為

式中，ti時刻的壓力增量ΔPi為

ti時刻的氣相泄放質量mv,i的差分關系(分亞聲速和超聲速兩種流態)為

式(10)～(12)的差分模型中，rb、χ和ρ可分別依照式(2)、式(5)和式(8)計算，而泄壓面積Ai可以簡化為分段函數

2.2 差分方程的求解步驟

確定參數后，在確定的初始條件下便可通過差分方程模型求解壓力變化過程和泄壓條件下的最大超壓(P-P0)max，模型中的壓力均為絕對壓力。其迭代求解算法如圖1所示。

圖1 差分方程求解流程框圖

3 模型的驗證與分析

3.1 時間差分的敏感性分析

差分離散的求解需要進行時間差分的敏感性分析，在本模型中即是針對時間差分的大小作敏感分析，以確定合適的Δt。現給定物質為甲烷，場景參數如表1所示，令Δt=0.01,0.005,0.002 5,0.001,0.000 5 s進行求解，其結果如圖2所示。由圖2可以看出：Δt的取值不同，模型求解的結果有一定偏差，當Δt<0.001 s之后，模型的求解結果已基本不受Δt取值的影響而趨于穩定。因此下文計算中，以Δt=0.001 s進行迭代求解可以得到穩定結果。

3.2 模型的可靠性分析

為驗證模型的可靠性，以Bartknecht在1 m3且有0.36 m2泄壓口的正方體空間內進行甲烷爆炸試驗條件進行計算對比[11]。試驗的實際條件參數如表2所示。為作區分，下文中用到的Pred指的是Bartknecht試驗下測得的超壓，而模型所計算的超壓均記為P-P0。

表1 敏感性分析場景參數表

圖2 求解的時間差分敏感性分析結果

由于泄壓構件的形式未知，因此假定Cd=0.90，且認為當空間內部壓力達到泄壓構件靜態觸發壓力的瞬間，泄壓面積完全敞開至最大且不隨壓力的下降而減小；初始的火焰燃燒半徑為0.001 m(rb,0/a=0.001 61→0)。計算與試驗結果對比情況如表3所示，各工況下的壓力變化過程如圖3所示。

表2 Bartknecht甲烷爆炸試驗條件參數匯總表

表3 模型-試驗最大超壓結果對比

由圖3可見，無論Pstat大小如何，空間內部超壓最初均按照相同的增速增長模式變化，這反映了在泄壓構件未開啟之前封閉空間的壓力增長模式。當P-P0大于Pstat-P0時，曲線出現拐點且此處的導數不連續，這是由于此時Av的突變引起的，類似爆破片等形式泄壓構件的動作模式；對于泄爆門等開口面積會隨著壓差大小變化的構件，超壓曲線則不會有此突躍。由表3的對比結果，可說明該模型存在一定的適用范圍：對于1 m3且泄壓面積為0.36 m2空間內的甲烷爆燃，當Pstat-P0處于20 kPa以下時，模型計算相對誤差較大，處于不可接受的范圍；而當Pstat-P0達到一定程度(150 kPa以上)，則相對誤差又會有所增大。因此，有必要進行模型參數的敏感性分析，以確定影響該模型適用條件的相關因素。

圖3 Bartknecht甲烷爆炸試驗的模型求解結果

3.3 參數的敏感性及誤差分析

3.3.1 壓力增長方式分析

現就3.2分析存在的問題，以表2的空間及物質條件，以2為等比梯度，分別求解在Pstat-P0從2.5～320 kPa下甲烷泄壓爆炸所能產生的最大超壓，對應的壓力變化過程如圖4所示，結果如表4所示。由圖4可得：超壓主要是取決于泄壓后的過程。在Pstat-P0≤20 kPa時，超壓的變化過程對泄壓動作已不敏感，在泄壓口開啟之后壓力繼續增長的時間較長且增長壓力相對接近；隨著Pstat的增大，壓力開始衰減的時間與泄壓口動作的時間越來越接近，達到Pstat-P0≥320 kPa時，在泄壓動作之后壓力隨即衰減而無法繼續增長。此二因素是導致模型在Pstat過高或過低的工況下誤差增大的主要原因。

圖4 (Pstat-P0)系列等比梯度模型求解結果表4 等比梯度工況計算結果

計算參量等比梯度工況Pstat-P0/kPa2.5510204080160320Pmax-P0/kPa29.6532.8941.6852.6971.64107.48168324泄壓動作時間/s0.0860.1090.1460.1770.2130.2660.3250.396衰減時間/s0.420.3990.3750.3660.3650.360.3570.397fvent,max0.9998880.9998790.9998210.9997790.9995970.9993990.9988960.997959fvent,min0.9953720.9950330.9941780.9931290.9914810.9890490.9864160.982043

現假設某個時刻ti存在Pi=(+)fvent>Pi-1，即(1-fvent)/fvent<ΔPi/Pi-1。由表4可得，fvent的值均在0.99左右，因此可認為不等式左側為一常數，壓力的增減主要由不等式右側的增壓比決定。該增壓比與(P/P0)γ、{1-(P/P0)γ[(Pm-P)/(Pm-P0)]}2/3成正比，即與P成負相關關系。因此在其他參數為常數的條件下，只能依靠P增長到一定程度后才能使其壓力衰減。因此，在其他參數保證可靠的前提下，該等比梯度系列計算結果是合理的。

3.3.2 湍流因子的敏感性及誤差機理分析

相關研究表明，KG并不是一個常數，其大小會受到湍流作用的影響，在不同條件下的數值應不同。而本文模型使用的湍流因子是通過恒定的KG倒推得到，本身具有局限性；同時，模型中的其他參數也直接或間接與χ相關，故χ是模型誤差的主要來源。現就表2工況1和工況7，以0.2的等差梯度研究χ=0.2～1.2的范圍內超壓的變化過程，得到的結果如表5、6所示，超壓變化過程如圖5、6所示。

表5 Pstat-P0=10 kPa工況下χ的敏感性分析結果

表6 Pstat-P0=200 kPa工況下χ的敏感性分析結果

圖5 Pstat-P0=10 kPa工況下χ的敏感性分析

圖6 Pstat-P0=200 kPa工況下χ的敏感性分析

參數χ體現的是爆炸中湍流程度的大小。從分析結果可以得到規律：當χ相對偏小時，泄壓口的動作時間和超壓開始衰減的時間將更接近，并可能在泄壓構件動作時壓力即開始衰減；當χ相對偏大時，泄壓口的動作時間和超壓開始衰減的時間將延長，造成泄壓后的超壓增長。對比結果：對于工況1，χ=0.6吻合程度較高；對于工況7,χ=1.2仍相對較小。說明Pstat-P0相對較低時，對于實際湍動程度給定的χ過大；而Pstat-P0相對較高時，χ又無法達到實際湍動水平。對該現象的解釋為：壓力的升高會導致傳播速度加快，引起湍動程度增大，所以較大的Pstat-P0會引起泄壓前維持較高內部壓力，進而對應的χ也應增大。NFPA 68中考慮泄壓引起的湍流程度時，所給雷諾數的計算模型也考慮到了不同壓力下泄放速度不同導致湍動程度不同的關系，并給出Rev∝P1/2的關系，這也印證了該誤差分析的合理性。

4 結論

研究可得到如下結論：(1)基于封閉空間的壓力增長理論和泄放質量比例的壓力衰減假設，建立了泄壓條件下氣體爆炸超壓變化的微分方程模型，并通過固定時間差分離散將其轉化為差分方程模型，給出了求解步驟。(2)模型進行代數求解時，時間差分的大小應控制在0.001 s以下以保證求解結果的穩定性。(3)恒定的湍流因子是導致模型誤差的主要原因，該模型進行泄壓條件下封閉空間的超壓計算時，應當考慮通過燃爆指數折算的湍流因子的適用范圍。體積為1 m3且泄壓面積為0.36 m2空間內的甲烷爆炸，在Pstat-P0處于20～200 kPa之間的誤差處于可接受水平。

值得一提的是，湍流程度受到諸多因素的影響，除上文分析的Pstat以外，不同的容器體積、泄壓面積、物質種類、初始壓力等條件下得到的χ也有所不同，因此模型適用范圍也不相同。針對此問題，可以通過后續研究湍流因子在泄壓條件下氣體爆炸過程中的變化規律進行補充。

[1] 中華人民共和國公安部.建筑設計防火規范:GB 50016—2014[S].北京:中國計劃出版社,2015.

[2] Standard on explosion protection by deflagration venting:NFPA 68—2013[S].

[3] Gas explosion venting protective system:BS EN 14994—2007[S].

[4] SWIFT I.Design of deflagration protection systems[J].Journal of Loss Prevention in the Process Industries,1988(1):5-15.

[5] BARTKNECHT W.Explosions-schutz:grundlagen and anwendung[M].Springer-Verlag,Berlin,1993.

[6] ZALOSH R.SFPE handbook of fire protection engineering [M].5th ed. Springer New York Heidelberg Dordrecht London:Society of Fire Protection Engineers,2016.

[7] BRADLEY M,MITCHESON A.The venting of gaseous explosions in spherical vessels[J].Combustion and Flame,1978,32:221-236.

[8] MOLKOV V,KOROLCHENKO A,ALEANDROV S.Venting of deflagrations in buildings and equipment:universal correlation[R].Russian Research Institute for Fire Protection.

[9] HERNANDEZ F,ABDEL-JAWAD M,HAO H.Simplified multiple equations’ inverse problem of vented vessels subjected to internal gas explosion[J].Journal of Loss Prevention in the Process Industries,2015,35:65-79.

[10] CHYZ Y T,MACKIEWICZ M.Simplified function of indoor gas explosion in residential building[J].Fire Safety Journal,2017,87:1-9.

[11] KASMAN RAFIZIANA M D.Vented gas explosion[M].School of Process,Environmental and Materials Engineering,University of Leeds,2008.

(責任編輯馬龍)

GasExplosionOverpressureCalculationModelunderDeflagrationVenting

YUAN Dancan, MA Bin, DU Kunze

(TeamofGraduateStudents,TheArmedPoliceAcademy,Langfang,HebeiProvince065000,China)

Based on the theory of pressure increase in enclosed space and the assumed pressure decrease related to the discharge mass ratio, a series of differential equation models are established to interpret the process of gas explosion overpressure under deflagration venting. Through difference discrete within fixed time, the differential equation models are transformed into the difference equation models with solution method proposed. The maximum time element of 0.001 s is suggested through the sensibility analysis of the time difference. Based on methane explosion experiments in a 1m3enclosure with a 0.36 m2vent conducted by Bartknecht, a series of simulations are computed for comparison. Additionally, the scope of application, the causes of error and the mechanism of turbulent factor are analyzed. The results show that when the relative static activation pressure of vent is between 20 to 200 kPa, Bartknecht’s experiments are reliable. When applied to other conditions of deflagration venting, the application range of turbulence factor calculated by the combustion index should be considered.

gas explosion; deflagration venting; overpressure model; turbulent factor

2017-05-21

袁丹燦(1993— )，男，廣東汕頭人，在讀碩士研究生； 馬斌(1994— )，男，甘肅蘭州人，在讀碩士研究生； 杜坤澤(1994— )，男，四川遂寧人，在讀碩士研究生。

D631.6

A

1008-2077(2017)10-0014-06