陸胤臣，陶 剛，張禮敬

（南京工業(yè)大學(xué)安全科學(xué)與工程學(xué)院，江蘇南京210009）

球形容器內(nèi)甲烷－空氣爆炸特性分析與理論計(jì)算＊

陸胤臣，陶 剛，張禮敬

（南京工業(yè)大學(xué)安全科學(xué)與工程學(xué)院，江蘇南京210009）

為研究甲烷－空氣混合物在密閉球形容器內(nèi)的爆炸特性，首先利用化學(xué)平衡計(jì)算軟件確定合適的燃燒產(chǎn)物與化學(xué)平衡溫度，估算甲烷－空氣混合物的最大爆炸壓力。然后基于火焰增長(zhǎng)模型，用MATLAB編輯循環(huán)語(yǔ)句程序，計(jì)算了甲烷－空氣混合氣體爆炸的壓力時(shí)程曲線，通過(guò)與實(shí)驗(yàn)數(shù)據(jù)對(duì)比，驗(yàn)證了化學(xué)平衡軟件計(jì)算方法與火焰增長(zhǎng)模型的可行性，并分析了誤差的產(chǎn)生原因。進(jìn)而利用火焰增長(zhǎng)模型推出的經(jīng)驗(yàn)公式計(jì)算爆燃指數(shù)，發(fā)現(xiàn)在當(dāng)量比附近與實(shí)驗(yàn)結(jié)果擬合程度較好。

最大爆炸壓力；爆燃指數(shù)；化學(xué)平衡計(jì)算；火焰增長(zhǎng)模型

天然氣作為目前廣泛使用的燃料，其主要成分甲烷與空氣的混合氣體在密閉球形容器內(nèi)的爆炸特性，如最大爆炸壓力、爆燃指數(shù)等是非常重要的物質(zhì)參數(shù)，可以用于風(fēng)險(xiǎn)分析與泄放設(shè)計(jì)。采用實(shí)驗(yàn)方法獲得這些數(shù)據(jù)需要昂貴的實(shí)驗(yàn)器材、訓(xùn)練有素的操作員和大量的時(shí)間，因此，進(jìn)行相關(guān)的理論和計(jì)算研究工作具有重要意義。

國(guó)內(nèi)外研究者已經(jīng)在爆炸特性方面進(jìn)行了較多的研究工作。D．Bradley等［1］建立了同時(shí)包含已燃區(qū)域和未燃區(qū)域的兩區(qū)域模型，為后來(lái)的理論和實(shí)驗(yàn)工作奠定了基礎(chǔ)；A．E．Dahoe等［2］研究了密閉容器中的粉塵爆炸，推導(dǎo)了粉塵爆燃指數(shù)Kst的公式，對(duì)可燃性氣體爆燃指數(shù)KG的推導(dǎo)有很大參考意義；C．V．Mashuga等［3］利用絕熱火焰溫度對(duì)爆炸極限進(jìn)行預(yù)測(cè)，為爆炸極限和最大爆炸壓力的研究提供了新思路；E．V．D．Bulck［4］在兩區(qū)域模型基礎(chǔ)上推導(dǎo)出了爆燃指數(shù)的計(jì)算公式，該公式形式簡(jiǎn)單，并明確給出爆燃指數(shù)與火焰燃燒速度聯(lián)系；Y．D．Jo等［56］建立了火焰增長(zhǎng)模型，并將其用于H2和空氣混合氣體壓力時(shí)程曲線的計(jì)算，發(fā)現(xiàn)計(jì)算結(jié)果與實(shí)驗(yàn)符合較好，進(jìn)而推導(dǎo)出爆燃指數(shù)的表達(dá)式；G．Tao等［7］利用化學(xué)平衡計(jì)算方法對(duì)H2／O2／N2混合氣體的爆炸壓力進(jìn)行了估算，通過(guò)與實(shí)驗(yàn)結(jié)果的對(duì)比發(fā)現(xiàn)，該結(jié)果具有較好的預(yù)測(cè)效果。

本文中，首先使用化學(xué)平衡方法計(jì)算得到最大爆炸壓力，然后基于火焰增長(zhǎng)模型，根據(jù)收集的參數(shù)數(shù)據(jù)，用Tablecurve2D擬合經(jīng)驗(yàn)公式，并利用MATLAB編寫程序，通過(guò)理論計(jì)算得到甲烷－空氣混合氣體的壓力時(shí)程曲線和爆燃指數(shù)，以期為最大爆炸壓力及壓力時(shí)程曲線的快速預(yù)測(cè)提供簡(jiǎn)便方法。

1 化學(xué)平衡計(jì)算

化學(xué)平衡計(jì)算方法是近年來(lái)出現(xiàn)的預(yù)測(cè)爆炸極限和最大爆炸壓力的方法，在氫氣研究中已經(jīng)取得較好的效果。爆炸平衡參數(shù)可以利用STANJAN、GASEQ等化學(xué)平衡計(jì)算軟件得到，假設(shè)燃燒滿足絕熱和固定體積條件，通過(guò)計(jì)算可得平衡態(tài)時(shí)的壓力、溫度和組分。本文中使用D．S．Dandy教授制作的在線軟件進(jìn)行化學(xué)平衡計(jì)算［8］，設(shè)初始?jí)毫?01 325Pa，初始溫度為298K。在該條件下，估算得到的平衡溫度為1 450K，平衡壓力為689．48kPa，接近實(shí)際的最大爆炸壓力。J．G．Du等［9］指出1 450K是可燃性混合物的合適限值，即絕熱火焰溫度高于該值時(shí)可以認(rèn)為是可燃的，低于該值則認(rèn)為是不可燃的，與我們的估算相符。在一定范圍內(nèi)，平衡壓力估算值的改變不會(huì)對(duì)計(jì)算結(jié)果造成明顯的影響。

計(jì)算最大的難點(diǎn)在于產(chǎn)物的確定，通過(guò)反復(fù)實(shí)驗(yàn)消除振蕩影響后發(fā)現(xiàn)，反應(yīng)后12種產(chǎn)物包括CH4、N2、O2、C、H、O、N、OH、CO、CO2、H2O、H2，其中C和H2對(duì)計(jì)算結(jié)果有較大的影響，必須列入產(chǎn)物中，而CH2、CH2OH、CH3O、H2O2、NH等物質(zhì)的含量極低，只有10－9量級(jí)，將其加入到產(chǎn)物中或去除對(duì)結(jié)果影響較小，甚至沒有影響?？諝饩€上化學(xué)平衡計(jì)算結(jié)果與實(shí)驗(yàn)值的比較如圖1所示，其中“p”為壓力，“α”為甲烷和空氣物質(zhì)的量之比（以下稱為空氣燃料比）。從圖1中可以看出：實(shí)驗(yàn)結(jié)果在當(dāng)量比附近（相當(dāng)于在α＝1．0附近取值）與計(jì)算結(jié)果符合較好；而在上、下限附近，兩者差別較大。這主要是因?yàn)楫?dāng)量比附近化學(xué)反應(yīng)速度快，火焰?zhèn)鞑パ杆伲磻?yīng)完全，與估計(jì)的產(chǎn)物相符程度較高；而在上、下限部分因?yàn)槿紵怀浞?，燃料沒有燃完，反應(yīng)速率慢，化學(xué)產(chǎn)物成分復(fù)雜，不易估計(jì)。另一重要的原因是重力對(duì)火焰?zhèn)鞑ズ蛪毫Φ挠绊?，在不考慮重力或微重力情況下，對(duì)所有可燃預(yù)混氣體，火焰都能基本保持球形對(duì)稱向外傳播，在上、下限附近火焰微弱，火焰形狀容易受到影響。在實(shí)際重力條件下，該情況會(huì)更加明顯，當(dāng)浮力引起的氣流上升速度大于燃燒速度時(shí)，火焰只能向上傳播，造成容器下半部分的燃料無(wú)法燃燒，造成壓力比理想狀態(tài)低很多，成分復(fù)雜。K．L．Cashdollar等［10］也得到了相似的結(jié)果。

圖1 化學(xué)平衡計(jì)算所得爆炸壓力值與實(shí)驗(yàn)的對(duì)比Fig．1 Comparison of chemical equilibrium calculation results with experimental data

2 參數(shù)確定及壓力時(shí)程曲線的理論計(jì)算

2．1 火焰增長(zhǎng)模型

火焰增長(zhǎng)模型是一種描述球形容器內(nèi)的爆炸過(guò)程的理論。該模型假設(shè)在中心點(diǎn)火后，火焰前沿以球形向外擴(kuò)展，并且以一個(gè)很小的微元步長(zhǎng)一層一層向外反應(yīng)，直至火焰前沿達(dá)到容器壁。與之前模型不同的是，該模型區(qū)分已燃部分和未燃部分的絕熱指數(shù)，火焰增長(zhǎng)情況更加接近實(shí)際燃燒過(guò)程，預(yù)測(cè)的壓力值也更加準(zhǔn)確?；鹧嬖鲩L(zhǎng)模型的假設(shè)如下：（1）火焰前沿為球形；（2）氣體壓力均一；（3）氣體浮力和火焰厚度可以忽略。基于以上假設(shè)，火焰增長(zhǎng)模型中壓力、燃燒爆炸程度、任意時(shí)刻半徑可由以下方程得到，即：

2．2 層流燃燒速率的確定

火焰增長(zhǎng)模型中，層流火焰速度Su的定義為平坦的、無(wú)拉伸的、絕熱的、預(yù)混的火焰前鋒傳播到未燃?xì)怏w時(shí)的火焰速度。層流火焰速度是燃燒混合物的重要參數(shù)，決定著燃燒程度和實(shí)際燃燒過(guò)程中火焰的穩(wěn)定性。層流火焰速度是化學(xué)機(jī)理研究的重要對(duì)象，許多科學(xué)家對(duì)甲烷／空氣混合物的層流燃燒速率進(jìn)行過(guò)測(cè)量。本文中對(duì)有代表性的數(shù)據(jù)進(jìn)行分析［1114］，并通過(guò)數(shù)值擬合得出層流火焰速度與空氣燃料比α之間的關(guān)系。

為使擬合公式具有工程應(yīng)用價(jià)值，本文使用三次多項(xiàng)式進(jìn)行擬合，所得表達(dá)式為y＝A＋B1x＋B2x2＋B3x3，參數(shù)的取值如下：A＝－183．12，B1＝448．52，B2＝－256．91，B3＝27．58。決定系數(shù)R2可以達(dá)到0．96954，即擬合效果較好。

2．3 絕熱指數(shù)的確定

絕熱指數(shù)是求解氣體膨脹程度的重要參數(shù)，在火焰增長(zhǎng)模型中，該參數(shù)在爆炸程度和壓力時(shí)程曲線計(jì)算中有很大作用。其他模型將已燃部分和未燃部分的絕熱指數(shù)取為相等的值，以簡(jiǎn)化計(jì)算，與實(shí)際情況有較大的不同，會(huì)造成計(jì)算上的誤差。E．V．D．Bulck［15］通過(guò)實(shí)驗(yàn)方法得到了甲烷已燃和未燃部分絕熱指數(shù)數(shù)據(jù)?；谖墨I(xiàn)［15］數(shù)據(jù)，本文利用Tablecurve 2D軟件進(jìn)行數(shù)值擬合，并從Tablecurve2D上千個(gè)擬合多項(xiàng)式中挑選出形式簡(jiǎn)單且擬合效果較好的結(jié)果，如圖3和圖4所示。從圖3和圖4可以看出；γu和γb擬合曲線的決定系數(shù)均非常接近1，說(shuō)明擬合效果很好。

圖2 層流燃燒速率與空氣燃料比擬合曲線Fig．2 Fitted curve of laminar flame speed vs．a(chǎn)ir fuel ratio

圖3 未燃部分絕熱指數(shù)擬合曲線Fig．3 Fitted curve of specific heat ratio for unburned gas

圖4 已燃部分絕熱指數(shù)擬合曲線Fig．4 Fitted curve of specific heat ratio for unburned gas

2．4 理論計(jì)算與比較

取空氣燃料比為1．055的情況進(jìn)行理論計(jì)算，根據(jù)擬合公式計(jì)算可得：已燃部分的絕熱指數(shù)約為1．062，未燃部分的絕熱指數(shù)為1．374，層流火焰速度根據(jù)擬合公式計(jì)算為36．51cm／s；利用化學(xué)平衡計(jì)算估計(jì)可得：絕熱狀態(tài)下最大壓力為912．87kPa，對(duì)應(yīng)的量綱一壓力為9．01；實(shí)驗(yàn)時(shí)每隔0．333ms記錄一個(gè)數(shù)據(jù)點(diǎn)，因此理論計(jì)算時(shí)間間隔取為0．333ms，對(duì)應(yīng)的量綱一時(shí)間為0．007 4；實(shí)驗(yàn)過(guò)程中點(diǎn)火方式為中心電火花點(diǎn)火，可以認(rèn)為點(diǎn)火區(qū)域集中于一個(gè)半徑很小的球形中，所以理論計(jì)算時(shí)點(diǎn)火半徑取為1 cm，對(duì)應(yīng)的量綱一半徑為0．06；甲烷的指數(shù)參數(shù)根據(jù)文獻(xiàn)［5－6］中數(shù)據(jù)取為0．1。用MATLAB進(jìn)行編程，所有的初始值和需要提前計(jì)算的數(shù)據(jù)均已給出，將其代入式（1）～式（3）進(jìn)行循環(huán)運(yùn)算，直到爆炸程度達(dá)到單位1。

圖5給出了基于火焰增長(zhǎng)模型計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)的比較。由圖5可知，當(dāng)爆炸程度達(dá)到1時(shí)，最大壓力為902．73kPa，與用化學(xué)平衡計(jì)算得到的爆炸壓力912．87kPa基本一致，也出現(xiàn)略高于實(shí)驗(yàn)值的情況，但更為精確。從圖3還可以看出，當(dāng)爆炸程度達(dá)到1時(shí)，量綱一半徑已經(jīng)略微超過(guò)1。這主要是因?yàn)樵谟肕ATLAB進(jìn)行迭代計(jì)算時(shí)必須給定一個(gè)很小的量綱一半徑初始值，而爆炸程度初始值為零，雖然初始值較小，但經(jīng)過(guò)接近300步的迭代和求根計(jì)算累積的誤差，造成量綱一半徑和爆炸程度不可能同時(shí)達(dá)到單位1；另一方面，火焰前沿接觸到容器壁時(shí)會(huì)出現(xiàn)拉伸、反射、壓縮的現(xiàn)象，會(huì)造成壓力進(jìn)一步的升高。總體來(lái)說(shuō)，甲烷－空氣混合物在空氣燃料比為1．055時(shí)，理論計(jì)算得到的壓力時(shí)程曲線與實(shí)驗(yàn)曲線的擬合效果較好，但略高于實(shí)驗(yàn)值。這主要是因?yàn)槔碚撚?jì)算假設(shè)是理想狀態(tài)，不考慮熱量散失、浮力等造成的。與Y．D．Jo等［56］所得結(jié)果相比，本文采用化學(xué)平衡計(jì)算預(yù)測(cè)最大壓力，假設(shè)平衡溫度為1 450K，參與反應(yīng)過(guò)程元素假設(shè)與其并不完全相同，部分元素對(duì)結(jié)果幾乎沒有影響。對(duì)于層流燃燒速率、已燃部分絕熱指數(shù)、未燃部分絕熱指數(shù)，本文總結(jié)了較多數(shù)據(jù)并進(jìn)行了擬合，有更好的工程應(yīng)用價(jià)值。此外，通過(guò)循環(huán)計(jì)算也能得到爆炸程度與時(shí)間的關(guān)系、爆炸壓力上升速度與時(shí)間的關(guān)系。

圖5 火焰增長(zhǎng)模型計(jì)算結(jié)果與實(shí)驗(yàn)值的比較Fig．5 Comparison of theoretical results by flame growth model with experimental data

3 爆燃指數(shù)的計(jì)算

應(yīng)用火焰增長(zhǎng)模型還可計(jì)算爆燃指數(shù)。爆燃指數(shù)是估計(jì)泄放面積、預(yù)防事故發(fā)生方面的重要參數(shù)，火焰增長(zhǎng)模型推導(dǎo)出的爆燃指數(shù)計(jì)算公式如下：

式（4）假設(shè)火焰前沿接觸爆炸容器球形表面時(shí)的狀態(tài)已知，未考慮爆炸結(jié)束段火焰表面積的變化。由于三次方定律，爆燃指數(shù)與體積無(wú)關(guān)，但事實(shí)上爆炸釜體積的變化會(huì)對(duì)爆燃指數(shù)產(chǎn)生影響。當(dāng)火焰前沿厚度不可忽略時(shí)，假設(shè)就會(huì)失效。

表1給出了不同爆炸程度條件下球形容器的爆燃指數(shù)，其中KG，exp為爆燃指數(shù)的實(shí)驗(yàn)值。由表1可知：當(dāng)ξ＝1時(shí)，利用式（4）計(jì)算得到的爆燃指數(shù)KG＝8．240 91MPa·m·s－1，比實(shí)驗(yàn)值高出較多；當(dāng)爆炸程度為0．95和0．90時(shí)，計(jì)算結(jié)果與實(shí)驗(yàn)值更加接近。上述結(jié)果表明，一部分的火焰前沿在爆炸完成（即爆炸程度ξ＝1）前就已經(jīng)接觸到了球形爆炸容器的表面，造成火焰表面積有很大程度的變化，對(duì)爆炸壓力上升速度有較大影響，從而影響爆燃指數(shù)值。因此，爆炸程度ξ＝0．90時(shí)，爆燃指數(shù)的計(jì)算值與實(shí)驗(yàn)值符合更好。

表1 球形容器內(nèi)甲烷－空氣混合氣體的爆燃指數(shù)Table 1 Deflagration index of methane－air mixture in spherical vessel

4 結(jié) 論

（1）利用化學(xué)平衡計(jì)算方法計(jì)算甲烷－空氣最大爆炸壓力時(shí)，在當(dāng)量比附近有較好的預(yù)測(cè)結(jié)果，在上、下限附近，由于燃燒機(jī)理的復(fù)雜性、燃燒速度較低，造成預(yù)測(cè)結(jié)果不理想。

（2）在當(dāng)量比附近，應(yīng)用火焰增長(zhǎng)模型得到的壓力時(shí)程曲線與實(shí)驗(yàn)數(shù)據(jù)符合性較好。

（3）應(yīng)用火焰增長(zhǎng)模型計(jì)算得到的爆燃指數(shù)KG與實(shí)驗(yàn)值符合性較好，并且使用爆炸程度為0．9時(shí)的壓力計(jì)算得到的爆燃指數(shù)更加準(zhǔn)確，說(shuō)明在釜中燃料燃燒結(jié)束前，爆炸壓力上升速度已達(dá)最大值。

感謝江蘇研究生科研創(chuàng)新計(jì)劃項(xiàng)目及Michigan Technological University的資助，同時(shí)感謝Daniel A．Crowl教授的幫助。

［1］ Bradley D，Mitcheson A．Mathematical solutions for explosions in spherical vessels［J］．Combustion and Flame，1976，26（2）：201－217．

［2］ Dahoe A E，Zevenbergen J F，Lemkowitz S M，et al．Dust explosion in spherical vessels：The role of flame thickness in the validity of the cube－root law［J］．Journal of Loss Prevention in the Process Industries，1996，18（9）：33－44．

［3］ Mashuga C V，Crowl D A．Flammability zone prediction using calculated adiabatic flame temperatures［J］．Process Safety Progress，1999，18（3）：127－134．

［4］ Bulck E V D．Closed algebraic expressions for the adiabatic limit value of the explosion constant in closed volume combustion［J］．Journal of Loss Prevention in the Process Industries，2005，18（1）：35－42．

［5］ Jo Y D，Crowl D A．Flame growth model for confined gas explosion［J］．Process Safety Progress，2009，28（2）：141－146．

［6］ Jo Y D，Crowl D A．Explosion characteristics of hydrogen－air mixtures in a spherical vessel［J］．Process Safety Progress，2010，29（3）：216－223．

［7］ Tao G，Crowl D A．Comparison of the maximum gas combustion pressure of hydrogen／oxygen／nitrogen between chemical equilibrium calculations and experimental data［J］．Procedia Engineering，2013，62：786－790．

［8］ Dandy D S．Bioanalytical microfluidics program［EB／OL］．［2015－12－28］．http：∥navier．engr．colostate．edu／tools／equil．html．

［9］ Du J G，Ma H H，Qu Z W，et al．Prediction of methanes flammability using chemical equilibrium［J］．Process Safety Progress，2015，34（1）：31－35．

［10］ Cashdollar K L，Zlochower I A，Green G M，et al．Flammability of methane，propane，and hydrogen gases［J］．Journal of Loss Prevention in the Process Industries，2000，13（3）：327－340．

［11］ Van Maaren A，Thung D S，De Goey L P H．Measurement of flame temperature and adiabatic burning velocity of methane／air mixtures［J］．Combustion Science and Technology，1994，96（4／5／6）：327－344．

［12］ Dahoe A E，De Goey L P H．On the determination of the laminar burning velocity from closed vessel gas explosions［J］．Journal of Loss Prevention in the Process Industries，2003，16（6）：457－478．

［13］ Varea E，Modica V，Vandel A，et al．Measurement of laminar burning velocity and Markstein length relative to fresh gases using a new postprocessing procedure：Application to laminar spherical flames for methane，ethanol and isooctane／air mixtures［J］．Combustion and Flame，2012，159（2）：577－590．

［14］ Chen Z．On the accuracy of laminar flame speeds measured from outwardly propagating spherical flames：Methane／air at normal temperature and pressure［J］．Combustion and Flame，2015，162（6）：2442－2453．

［15］ Bulck E V D．Closed algebraic expressions for the adiabatic limit value of the explosion constant in closed volume combustion［J］．Journal of Loss Prevention in the Process Industries，2005，18（1）：35－42．

［16］ Benedetto A D，Cammarota F，Sarli V D，et al．Anomalous behavior during explosions of CH4in oxygen－enriched air［J］．Combustion and Flame，2011，158（11）：2214－2219．

Analysis and theoretical calculation of explosion characteristics of methane－air mixture in a spherical vessel

Lu Yinchen，Tao Gang，Zhang Lijing

（College of Safety Science and Engineering，Nanjing Tech University，Nanjing210009，Jiangsu，China）

To study the characteristics of the methane－air mixture exploding in a closed spherical container，we determined the appropriate combustion products and chemical equilibrium temperature using the chemical equilibrium calculation software，thereby predicting the maximum explosion pressure of the mixture．The MATLAB program based on the flame growth model was adopted to calculate the curve showing the relationship between the explosion pressure and time．The calculation processes were verified by the comparison of the obtained results with the experimental data，and the origin of the error was also identified．Further，it is found that the empirical formula of the deflagration index KGderived from the flame growth model is well fitted with the experimental date near the chemical equivalent line．

maximum pressure；deflagration index；chemical equilibrium calculation；flame growth model

O354．1國(guó)標(biāo)學(xué)科代碼：1303510

A

10．11883／1001－1455（2017）04－0773－06

（責(zé)任編輯 王玉鋒）

2015－12－28；

2016－03－30

江蘇省2015年度普通高校研究生科研創(chuàng)新計(jì)劃項(xiàng)目（KYLX15＿0774）

陸胤臣（1991－ ），男，碩士研究生；通信作者：陶 剛，taogang＠njtech．edu．cn。