徐維錚，吳衛(wèi)國(guó)，況 正

(1. 武漢理工大學(xué) 高性能艦船技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430063；2. 武漢理工大學(xué) 交通學(xué)院船舶、海洋與結(jié)構(gòu)工程系，湖北 武漢 430063)

0 引 言

工業(yè)可燃性氣云爆炸事故大多是由弱點(diǎn)火點(diǎn)燃可燃?xì)怏w引起的，其傳播形式大多是以亞音速傳播的爆燃波。爆燃波的傳播過(guò)程很復(fù)雜，受環(huán)境條件和物理因素的影響極大，當(dāng)氣體在傳播過(guò)程中由于障礙物阻礙發(fā)生湍流加強(qiáng)效應(yīng)或是在初始時(shí)刻由高能量直接起爆，其爆炸模式將由爆燃轉(zhuǎn)為爆轟或直接以爆轟模式發(fā)生爆炸。氣體發(fā)生爆轟以后將在空氣中形成強(qiáng)度較大的帶有負(fù)壓區(qū)的空氣沖擊波，對(duì)工作人員和周邊結(jié)構(gòu)設(shè)施將造成較大的損害。針對(duì)空中可燃?xì)庠?、云霧爆炸模擬，近年來(lái)徐勝利等學(xué)者進(jìn)行過(guò)相應(yīng)的數(shù)值研究[1 - 4]，他們將可燃?xì)庠啤⒃旗F等效的簡(jiǎn)化為一團(tuán)高壓氣體，而艙室內(nèi)部可燃?xì)庠票ㄑ芯枯^少。

爆炸屬于強(qiáng)間斷問(wèn)題，其數(shù)值模擬對(duì)激波捕捉格式提出較嚴(yán)格要求。目前，WENO格式作為一種典型的高精度激波捕捉格式，對(duì)流場(chǎng)內(nèi)的激波間斷具有較高的分辨率，適于求解包含激波、膨脹波以及接觸間斷等復(fù)雜結(jié)構(gòu)的流場(chǎng)。LIU等[5]于1994年提出了 WENO（weighted essentially non-oscillation scheme）格式，之后，SHU 等[6 - 7]發(fā)展了該格式并擴(kuò)展了其應(yīng)用。

考慮到WENO格式高精度、穩(wěn)定性較好的優(yōu)勢(shì)，基于FORTRAN平臺(tái)，采用三階WENO有限差分格式，自主開(kāi)發(fā)了密閉多艙室內(nèi)爆炸波二維數(shù)值計(jì)算程序。利用所開(kāi)發(fā)的程序開(kāi)展了多艙室內(nèi)柱形氣云爆炸波傳播過(guò)程及不同艙室內(nèi)部爆炸載荷特性研究。

1 控制方程

含激波流場(chǎng)采用可壓縮歐拉方程進(jìn)行描述，其具體形式如下：

在方程（1）的每個(gè)方向上均可以看成是一個(gè)雙曲守恒律方程：

例如，針對(duì)x方向，方程（5）的數(shù)值離散形式為：

2 數(shù)值方法

2.1 三階 WENO-JS 格式

三階WENO格式的數(shù)值離散和推導(dǎo)過(guò)程如下[6]：

利用上述2種二階通量的凸組合計(jì)算最終具有三階精度的數(shù)值通量，即

其中，本文參數(shù)ε取值為10-6。光滑因子的表達(dá)式如下：

2.2 時(shí)間項(xiàng)數(shù)值離散

采用三階TVD-RK法對(duì)歐拉方程（1）中的時(shí)間項(xiàng)進(jìn)行數(shù)值離散，具體的離散格式如下[6]：

3 試驗(yàn)驗(yàn)證

為了檢驗(yàn)所開(kāi)發(fā)數(shù)值程序的可靠性和正確性，選取文獻(xiàn)[8]開(kāi)展的平面激波繞方腔流動(dòng)試驗(yàn)進(jìn)行對(duì)比分析。圖1給出計(jì)算域的尺寸和網(wǎng)格劃分示意圖（本文計(jì)算中，網(wǎng)格尺寸為1 mm）。

在標(biāo)準(zhǔn)狀況下，平面激波馬赫數(shù)為1.3，根據(jù)沖擊波物理可以確定波后物理參數(shù)。圖2給出該算例初始條件。

圖3給出不同時(shí)刻數(shù)值模擬紋影圖與試驗(yàn)的對(duì)比，圖中不同時(shí)刻的波系結(jié)構(gòu)和臺(tái)階處繞射形成的漩渦結(jié)構(gòu)均證明了文中程序開(kāi)發(fā)的可靠性和正確性，為進(jìn)行多艙室內(nèi)爆炸模擬奠定了基礎(chǔ)。

圖1 計(jì)算區(qū)域尺寸及網(wǎng)格劃分Fig.1 Size of the computational zone

圖2 激波繞方腔流動(dòng)初始條件Fig.2 Initial condition of the planar shock wave over a square cavity

圖3 數(shù)值與文獻(xiàn)[8]中試驗(yàn)結(jié)果對(duì)比Fig.3 Comparisons between numerical and experimental results[8]

4 多艙室內(nèi)爆炸過(guò)程數(shù)值計(jì)算

本節(jié)主要采用驗(yàn)證的程序?qū)γ荛]多艙室內(nèi)平面柱形高壓氣云爆炸過(guò)程進(jìn)行數(shù)值模擬研究，主要探討爆炸波傳播路徑及不同艙室內(nèi)部爆炸載荷特性，為后期的結(jié)構(gòu)抗爆設(shè)計(jì)及毀傷評(píng)估提供一定的參考價(jià)值。

4.1 艙室尺寸及測(cè)點(diǎn)布置

多艙室有3個(gè)艙室組成，分別標(biāo)記為艙室1、艙室2、艙室3，艙室1與艙室2及艙室3之間通過(guò)方形連接導(dǎo)管相連通，每個(gè)艙室的具體尺寸參見(jiàn) 圖4。每個(gè)艙室的壁面中點(diǎn)處均設(shè)置2個(gè)壓力測(cè)點(diǎn)（分別編號(hào)為 NO1，NO2，NO3，NO4，NO5，NO6）對(duì)爆炸超壓載荷時(shí)間歷程進(jìn)行輸出（見(jiàn)圖4（a））。壁面邊界條件設(shè)置為剛性邊界（不考慮艙室結(jié)構(gòu)的變形），均勻網(wǎng)格尺寸取為1 mm（見(jiàn)圖4（b））。

圖4 多艙室布置及網(wǎng)格劃分Fig.4 Size of multi-cabin and mesh distribution

4.2 初始條件及工況設(shè)置

柱形均勻高壓氣云的具體參數(shù)：半徑為20 mm，密度為 1.99 kg/m3，壓力為 13.8×105Pa。為了探討爆炸發(fā)生位置對(duì)爆炸過(guò)程的影響規(guī)律，本文選取2種典型爆炸工況：工況1表示爆炸發(fā)生在艙室1內(nèi)部，如圖5（a）所示，工況2表示爆炸發(fā)生在艙室3內(nèi)部，如圖5（b）所示。均勻高壓氣云均設(shè)置在艙室中間（見(jiàn)圖5中的淺色區(qū)域），周圍深色區(qū)域表示空氣域，空氣密度為 1.225 kg/m3，壓力為 1.0×105Pa。

4.3 工況一爆炸過(guò)程分析

本節(jié)主要根據(jù)數(shù)值模擬結(jié)果分析工況1爆炸后爆炸波傳播過(guò)程及不同艙室爆炸載荷分布特性。

4.3.1 爆炸波傳播過(guò)程分析

為了探討工況1氣云爆炸后爆炸波在多艙室內(nèi)部的傳播過(guò)程，這里給出不同時(shí)刻爆炸壓力云圖和數(shù)值紋影圖（見(jiàn)圖6）。根據(jù)圖6可知工況1多艙室內(nèi)爆炸過(guò)程爆炸波的傳播路徑：爆炸波首先在艙室1中進(jìn)行自由膨脹運(yùn)動(dòng)（圖6（a））爆炸波初次到達(dá)艙室壁面時(shí)發(fā)生正規(guī)則反射，壁面處形成局部高壓區(qū)域（圖6（b）），之后爆炸波沿著連接導(dǎo)管進(jìn)入艙室2和艙室3內(nèi)部（圖6（c）），并緊接著在艙室2和艙室3上下壁面處進(jìn)行系列正規(guī)則斜反射和馬赫反射（圖6（d）和圖6（e）），當(dāng)爆炸波傳播到艙室2和艙室3的最右端壁面處時(shí)再次形成復(fù)雜的壁面反射波（圖6（f））。由于多艙室處于封閉空間，因此可以預(yù)知爆炸波經(jīng)過(guò)上述在不同艙室之間的傳播過(guò)程，爆炸波強(qiáng)度將不斷衰弱，爆炸氣體的動(dòng)能逐漸轉(zhuǎn)化為氣體的內(nèi)能，最終艙室內(nèi)部將趨近于一個(gè)準(zhǔn)靜態(tài)過(guò)程。

4.3.2 爆炸載荷特性分析

為了探討工況1氣云爆炸后不同艙室內(nèi)部爆炸載荷分布特性，這里給出不同艙室內(nèi)部壁面測(cè)點(diǎn)超壓時(shí)間歷程曲線和沖量歷程時(shí)間曲線。針對(duì)工況1，由于艙室2和艙室3上下對(duì)稱，因此，這里只需要討論艙室1和艙室2內(nèi)部壁面測(cè)點(diǎn)的爆炸載荷特性。圖7給出艙室1壁面測(cè)點(diǎn)NO1，NO2爆炸載荷時(shí)間歷程，從圖7可知艙室內(nèi)爆炸載荷主要呈現(xiàn)峰值不斷衰減的沖擊波載荷和準(zhǔn)靜態(tài)壓力。雖然測(cè)點(diǎn)NO1，NO2與爆源處于近似相同的距離，但是測(cè)點(diǎn)NO1超壓載荷及沖量明顯大于測(cè)點(diǎn)NO2，通過(guò)分析發(fā)現(xiàn)這主要是由于測(cè)點(diǎn)NO2靠近艙室2和艙室3，當(dāng)爆炸波通過(guò)連接導(dǎo)管傳播到艙室2和艙室3內(nèi)部后將明顯減弱測(cè)點(diǎn)NO2處的沖擊波強(qiáng)度。

圖5 初始條件及爆炸工況Fig.5 Initial condition and explosion cases

圖6 工況 1 不同時(shí)刻爆炸波傳播過(guò)程Fig.6 Blast wave propagation process at different time for explosion case one

圖7 艙室 1 壁面測(cè)點(diǎn) NO1，NO2 爆炸載荷時(shí)間歷程（工況1）Fig.7 Blast load time histories of gauging point NO1, NO2 in cabin 1 for explosion case one

圖8 給出艙室2壁面測(cè)點(diǎn)NO3，NO4爆炸載荷時(shí)間歷程，從圖8可知在爆炸初期，測(cè)點(diǎn)NO3超壓載荷及沖量稍微大于測(cè)點(diǎn)NO4，而在爆炸后期，測(cè)點(diǎn)NO3超壓載荷及沖量演變?yōu)樾∮跍y(cè)點(diǎn)NO4，通過(guò)分析發(fā)現(xiàn)這主要是由于爆炸初期測(cè)點(diǎn)NO3相較測(cè)點(diǎn)NO4更靠近爆源，因此其接觸到的初始沖擊波強(qiáng)度較大。而在爆炸后期由于測(cè)點(diǎn)NO4與連接導(dǎo)管正對(duì)，在后續(xù)的爆炸波傳播過(guò)程中，其將會(huì)受到爆炸波的連續(xù)直接沖擊，反而造成其超壓載荷、沖量時(shí)間歷程稍微高于測(cè)點(diǎn)NO3。

圖9給出艙室1，艙室2壁面測(cè)點(diǎn)NO1，NO2，NO3，NO4爆炸載荷時(shí)間歷程，從圖9可知艙室1內(nèi)部測(cè)點(diǎn)NO1，NO2爆炸載荷強(qiáng)度明顯強(qiáng)于艙室2內(nèi)部測(cè)點(diǎn)NO3，NO4，通過(guò)分析發(fā)現(xiàn)這主要是由于爆炸發(fā)生在艙室1內(nèi)部，因此其主要能量集中在艙室1內(nèi)部。同時(shí)根據(jù)圖9（b）中的沖量歷程曲線還得知，不同艙室內(nèi)部最終的準(zhǔn)靜態(tài)壓力相同，如果將沖量積分曲線開(kāi)始變?yōu)橹本€段的轉(zhuǎn)折點(diǎn)作為準(zhǔn)靜態(tài)壓力開(kāi)始形成的時(shí)間，則在該爆炸工況下不同艙室內(nèi)部測(cè)點(diǎn)到達(dá)準(zhǔn)靜態(tài)壓力的時(shí)間基本相同（近似為0.5 ms）。

圖8 艙室 2 壁面測(cè)點(diǎn) NO3，NO4 爆炸載荷時(shí)間歷程（工況1）Fig.8 Blast load time histories of gauging point NO3, NO4 in cabin 2 for explosion case one

圖9 艙室 1 及艙室 2 壁面測(cè)點(diǎn) NO1，NO2，NO3，NO4 爆炸載荷時(shí)間歷程（工況1）Fig.9 Blast load time histories of gauging point NO1, NO2, NO3,NO4 in cabin 1, 2 for explosion case one

圖10 工況 2 不同時(shí)刻爆炸波傳播過(guò)程Fig.10 Blast wave propagation process at different time for explosion case two

4.4 工況 2 爆炸過(guò)程分析

本節(jié)主要根據(jù)數(shù)值模擬結(jié)果分析工況2爆炸后爆炸波傳播過(guò)程及不同艙室爆炸載荷分布特性。

4.4.1 爆炸波傳播過(guò)程分析

為了探討工況2氣云爆炸后爆炸波在多艙室內(nèi)部的傳播過(guò)程，這里給出不同時(shí)刻爆炸壓力云圖和數(shù)值紋影圖（見(jiàn)圖10）。根據(jù)圖10可知工況2多艙室內(nèi)爆炸過(guò)程爆炸波的傳播路徑：爆炸波首先在長(zhǎng)方形艙室3中進(jìn)行自由膨脹運(yùn)動(dòng)，初次到達(dá)艙室上下壁面時(shí)發(fā)生正規(guī)則反射，壁面處形成局部高壓區(qū)域（圖10（a）），之后爆炸波沿著連接導(dǎo)管進(jìn)入艙室1內(nèi)部并以球形方式進(jìn)行擴(kuò)張（圖10（b）和圖10（c）），在擴(kuò)張過(guò)程中由于艙室1上壁面的限制，形成局部馬赫反射波（圖10（d）），緊接著，爆炸波傳播到艙室2內(nèi)部并形成復(fù)雜的壁面反射波（圖10（e）和圖10（f））。

圖11 艙室 1 壁面測(cè)點(diǎn) NO1，NO2 爆炸載荷時(shí)間歷程（工況2）Fig.11 Blast load time histories of gauging point NO1, NO2 in cabin 1 for explosion case two

圖12 艙室 2 壁面測(cè)點(diǎn) NO3，NO4 爆炸載荷時(shí)間歷程（工況2）Fig.12 Blast load time histories of gauging point NO3, NO4 in cabin 2 for explosion case two

4.4.2 爆炸載荷特性分析

為了探討工況2氣云爆炸后不同艙室內(nèi)部爆炸載荷分布特性，這里給出不同艙室內(nèi)部壁面測(cè)點(diǎn)超壓時(shí)間歷程曲線和沖量歷程時(shí)間曲線。圖11給出艙室1壁面測(cè)點(diǎn)NO1，NO2爆炸載荷時(shí)間歷程，從圖11可知測(cè)點(diǎn)NO1超壓載荷及沖量稍微大于測(cè)點(diǎn)NO2，通過(guò)分析發(fā)現(xiàn)這主要是由于測(cè)點(diǎn)NO1靠近艙室3的連接導(dǎo)管，當(dāng)爆炸波通過(guò)連接導(dǎo)管傳播到艙室1內(nèi)部時(shí)，測(cè)點(diǎn)NO1將遭受到更強(qiáng)的爆炸波沖擊。

圖12給出艙室2壁面測(cè)點(diǎn)NO3，NO4爆炸載荷時(shí)間歷程，從圖12可知測(cè)點(diǎn)NO3，NO4雖然超壓載荷時(shí)間歷程差異較大，但是沖量載荷歷程基本相同，通過(guò)分析發(fā)現(xiàn)這主要是由于測(cè)點(diǎn)NO3，NO4位于艙室2內(nèi)部，相較于艙室1中的爆源處于遠(yuǎn)場(chǎng)區(qū)域，當(dāng)爆炸波通過(guò)連接導(dǎo)管擴(kuò)張到艙室1再傳播到艙室2內(nèi)部后，由于艙室2容積也較艙室1偏小，傳入的爆炸波沖擊很快使得艙室2內(nèi)部達(dá)到準(zhǔn)靜態(tài)過(guò)程。

圖13給出艙室3壁面測(cè)點(diǎn)NO5，NO6爆炸載荷時(shí)間歷程，從圖13可知測(cè)點(diǎn)NO6初始超壓載荷明顯強(qiáng)于測(cè)點(diǎn)NO5，后期強(qiáng)度衰減速度快于測(cè)點(diǎn)NO5，這主要是由于艙室3為長(zhǎng)方形艙室，測(cè)點(diǎn)6距離爆源較近，因此初期沖擊波載荷強(qiáng)度大，由于測(cè)點(diǎn)6距離連接導(dǎo)管近，因此其后期超壓衰減速率更快。艙室3容積較小，雖然測(cè)點(diǎn)NO6距離連接導(dǎo)管較近，對(duì)該處的沖量歷程有一定的影響，但是測(cè)點(diǎn)6距離爆源較近，最終使得NO5，NO6沖量歷程曲線基本重合。

圖14給出艙室1～艙室3壁面測(cè)點(diǎn)NO1～NO6爆炸載荷時(shí)間歷程，從圖14可知艙室3內(nèi)部測(cè)點(diǎn)NO5，NO6爆炸載荷強(qiáng)度明顯強(qiáng)于艙室1內(nèi)部測(cè)點(diǎn)NO1，NO2，通過(guò)分析發(fā)現(xiàn)這主要是由于爆炸發(fā)生在艙室3內(nèi)部，因此其主要能量集中在艙室3內(nèi)部。艙室1內(nèi)部測(cè)點(diǎn)NO1，NO2爆炸載荷強(qiáng)度稍微強(qiáng)于艙室2內(nèi)部測(cè)點(diǎn)NO3，NO4，通過(guò)分析發(fā)現(xiàn)這主要是由于NO1，NO2相較測(cè)點(diǎn)NO3，NO4距離爆炸源更近一些，沖擊波強(qiáng)度隨著傳播距離的增大，強(qiáng)度逐漸衰弱。同時(shí)根據(jù)圖14（b）中的沖量歷程曲線還得知，盡管爆炸波達(dá)到各測(cè)點(diǎn)的時(shí)間存在差異性，但是該爆炸工況下不同艙室內(nèi)部測(cè)點(diǎn)到達(dá)準(zhǔn)靜態(tài)壓力的時(shí)間基本相同（近似為0.5 ms），這一結(jié)論與4.3.2小節(jié)中的分析結(jié)果相同。因此，通過(guò)這2個(gè)工況的爆炸分析，可以初步得出如下結(jié)論：對(duì)于特定的封閉艙室，準(zhǔn)靜態(tài)壓力形成時(shí)間主要與艙室容積和爆源能量相關(guān)。其實(shí)，從準(zhǔn)靜態(tài)壓力形成的物理過(guò)程來(lái)分析也很容易闡述這種關(guān)系，由于準(zhǔn)靜態(tài)壓力全場(chǎng)均勻，因此其在全場(chǎng)空間內(nèi)必然應(yīng)同時(shí)到達(dá)。

圖13 艙室 3 壁面測(cè)點(diǎn) NO5，NO6 爆炸載荷時(shí)間歷程（工況2）Fig.13 Blast load time histories of gauging point NO5, NO6 in cabin 3 for explosion case two

圖14 艙室 1～艙室 3 壁面測(cè)點(diǎn) NO1，NO2，NO3，NO4，NO5，NO6爆炸載荷時(shí)間歷程（工況2）Fig.14 Blast load time histories of gauging point NO1, NO2, NO3,NO4, NO5, NO6 in cabin 1,2,3 for explosion case two

圖15 工況1及工況2壁面測(cè)點(diǎn)NO1爆炸載荷時(shí)間歷程對(duì)比圖Fig.15 Comparisons of blast load time histories of gauging point NO1between explosion case one and two

4.5 工況 1 與工況 2 對(duì)比分析

為了進(jìn)一步分析爆源位置對(duì)多艙室內(nèi)部不同測(cè)點(diǎn)處爆炸載荷特性的影響規(guī)律，圖15～圖20分別給出2種工況下測(cè)點(diǎn)NO1～NO6的超壓載荷和沖量時(shí)間歷程曲線圖。

根據(jù)圖15和圖16可知，針對(duì)艙室1內(nèi)部測(cè)點(diǎn)NO1，NO2，工況1相較工況2具有更強(qiáng)的超壓載荷和更大的沖量，這主要是由于工況1中爆炸發(fā)生在艙室1中。

根據(jù)圖17和圖18可知，針對(duì)艙室2內(nèi)部測(cè)點(diǎn)NO3，NO4，2種工況超壓載荷歷程差異較大，但是沖量載荷歷程差別較小，這主要是由于艙室2中的測(cè)點(diǎn)相較2種工況中的爆源距離較大，屬于遠(yuǎn)場(chǎng)的范圍，沖量的歷程主要由準(zhǔn)靜態(tài)壓力決定。

根據(jù)圖19和圖20可知，針對(duì)艙室3內(nèi)部測(cè)點(diǎn)NO5，NO6，工況2相較工況1具有更強(qiáng)的超壓載荷和更大的沖量，這主要是由于工況2中爆炸發(fā)生在艙室3中。

圖16 工況1及工況2壁面測(cè)點(diǎn)NO2爆炸載荷時(shí)間歷程對(duì)比圖Fig.16 Comparisons of blast load time histories of gauging point NO2 between explosion case one and two

圖17 工況1及工況2壁面測(cè)點(diǎn)NO3爆炸載荷時(shí)間歷程對(duì)比圖Fig.17 Comparisons of blast load time histories of gauging point NO3 between explosion case one and two

圖18 工況1及工況2壁面測(cè)點(diǎn)NO4爆炸載荷時(shí)間歷程對(duì)比圖Fig.18 Comparisons of blast load time histories of gauging point NO4 between explosion case one and two

圖19 工況1及工況2壁面測(cè)點(diǎn)NO5爆炸載荷時(shí)間歷程對(duì)比圖Fig.19 Comparisons of blast load time histories of gauging point NO5 between explosion case one and two

圖20 工況1及工況2壁面測(cè)點(diǎn)NO5爆炸載荷時(shí)間歷程對(duì)比圖Fig.20 Comparisons of blast load time histories of gauging point NO5 between explosion case one and two

5 結(jié) 語(yǔ)

采用三階WENO有限差分格式，自主開(kāi)發(fā)了密閉多艙室內(nèi)爆炸波二維數(shù)值計(jì)算程序。利用所開(kāi)發(fā)的程序開(kāi)展了多艙室內(nèi)柱形氣云爆炸波傳播過(guò)程及不同艙室內(nèi)部爆炸載荷特性研究。通過(guò)本文的研究主要得到以下幾點(diǎn)結(jié)論：

1）多艙室內(nèi)爆炸過(guò)程中，爆炸所在艙室內(nèi)部爆炸載荷強(qiáng)度明顯強(qiáng)于相鄰艙室，尤其需要重點(diǎn)關(guān)注爆炸所在艙室通過(guò)連接導(dǎo)管對(duì)直接相鄰艙室爆炸載荷的影響；爆源位置對(duì)不同艙室內(nèi)部超壓載荷時(shí)間歷程影響較大，而對(duì)處于遠(yuǎn)場(chǎng)的艙室內(nèi)沖量時(shí)間歷程影響較小。

2）對(duì)于特定的封閉艙室，準(zhǔn)靜態(tài)壓力峰值及準(zhǔn)靜態(tài)壓力形成時(shí)間均全場(chǎng)均勻，準(zhǔn)靜態(tài)壓力形成時(shí)間主要與艙室容積和爆源能量相關(guān)。

3）本文研究可為工程抗爆結(jié)構(gòu)設(shè)計(jì)及爆炸毀傷評(píng)估提供一定的參考和指導(dǎo)。