徐維錚,吳衛(wèi)國,況 正
(1. 武漢理工大學 高性能艦船技術(shù)教育部重點實驗室,湖北 武漢 430063;2. 武漢理工大學 交通學院船舶、海洋與結(jié)構(gòu)工程系,湖北 武漢 430063)
工業(yè)可燃性氣云爆炸事故大多是由弱點火點燃可燃氣體引起的,其傳播形式大多是以亞音速傳播的爆燃波。爆燃波的傳播過程很復雜,受環(huán)境條件和物理因素的影響極大,當氣體在傳播過程中由于障礙物阻礙發(fā)生湍流加強效應或是在初始時刻由高能量直接起爆,其爆炸模式將由爆燃轉(zhuǎn)為爆轟或直接以爆轟模式發(fā)生爆炸。氣體發(fā)生爆轟以后將在空氣中形成強度較大的帶有負壓區(qū)的空氣沖擊波,對工作人員和周邊結(jié)構(gòu)設(shè)施將造成較大的損害。針對空中可燃氣云、云霧爆炸模擬,近年來徐勝利等學者進行過相應的數(shù)值研究[1 - 4],他們將可燃氣云、云霧等效的簡化為一團高壓氣體,而艙室內(nèi)部可燃氣云爆炸研究較少。
爆炸屬于強間斷問題,其數(shù)值模擬對激波捕捉格式提出較嚴格要求。目前,WENO格式作為一種典型的高精度激波捕捉格式,對流場內(nèi)的激波間斷具有較高的分辨率,適于求解包含激波、膨脹波以及接觸間斷等復雜結(jié)構(gòu)的流場。LIU等[5]于1994年提出了 WENO(weighted essentially non-oscillation scheme)格式,之后,SHU 等[6 - 7]發(fā)展了該格式并擴展了其應用。
考慮到WENO格式高精度、穩(wěn)定性較好的優(yōu)勢,基于FORTRAN平臺,采用三階WENO有限差分格式,自主開發(fā)了密閉多艙室內(nèi)爆炸波二維數(shù)值計算程序。利用所開發(fā)的程序開展了多艙室內(nèi)柱形氣云爆炸波傳播過程及不同艙室內(nèi)部爆炸載荷特性研究。
含激波流場采用可壓縮歐拉方程進行描述,其具體形式如下:
在方程(1)的每個方向上均可以看成是一個雙曲守恒律方程:
例如,針對x方向,方程(5)的數(shù)值離散形式為:
三階WENO格式的數(shù)值離散和推導過程如下[6]:
利用上述2種二階通量的凸組合計算最終具有三階精度的數(shù)值通量,即
其中,本文參數(shù)ε取值為10-6。光滑因子的表達式如下:
采用三階TVD-RK法對歐拉方程(1)中的時間項進行數(shù)值離散,具體的離散格式如下[6]:
為了檢驗所開發(fā)數(shù)值程序的可靠性和正確性,選取文獻[8]開展的平面激波繞方腔流動試驗進行對比分析。圖1給出計算域的尺寸和網(wǎng)格劃分示意圖(本文計算中,網(wǎng)格尺寸為1 mm)。
在標準狀況下,平面激波馬赫數(shù)為1.3,根據(jù)沖擊波物理可以確定波后物理參數(shù)。圖2給出該算例初始條件。
圖3給出不同時刻數(shù)值模擬紋影圖與試驗的對比,圖中不同時刻的波系結(jié)構(gòu)和臺階處繞射形成的漩渦結(jié)構(gòu)均證明了文中程序開發(fā)的可靠性和正確性,為進行多艙室內(nèi)爆炸模擬奠定了基礎(chǔ)。
圖1 計算區(qū)域尺寸及網(wǎng)格劃分Fig.1 Size of the computational zone
圖2 激波繞方腔流動初始條件Fig.2 Initial condition of the planar shock wave over a square cavity
圖3 數(shù)值與文獻[8]中試驗結(jié)果對比Fig.3 Comparisons between numerical and experimental results[8]
本節(jié)主要采用驗證的程序?qū)γ荛]多艙室內(nèi)平面柱形高壓氣云爆炸過程進行數(shù)值模擬研究,主要探討爆炸波傳播路徑及不同艙室內(nèi)部爆炸載荷特性,為后期的結(jié)構(gòu)抗爆設(shè)計及毀傷評估提供一定的參考價值。
多艙室有3個艙室組成,分別標記為艙室1、艙室2、艙室3,艙室1與艙室2及艙室3之間通過方形連接導管相連通,每個艙室的具體尺寸參見 圖4。每個艙室的壁面中點處均設(shè)置2個壓力測點(分別編號為 NO1,NO2,NO3,NO4,NO5,NO6)對爆炸超壓載荷時間歷程進行輸出(見圖4(a))。壁面邊界條件設(shè)置為剛性邊界(不考慮艙室結(jié)構(gòu)的變形),均勻網(wǎng)格尺寸取為1 mm(見圖4(b))。
圖4 多艙室布置及網(wǎng)格劃分Fig.4 Size of multi-cabin and mesh distribution
柱形均勻高壓氣云的具體參數(shù):半徑為20 mm,密度為 1.99 kg/m3,壓力為 13.8×105Pa。為了探討爆炸發(fā)生位置對爆炸過程的影響規(guī)律,本文選取2種典型爆炸工況:工況1表示爆炸發(fā)生在艙室1內(nèi)部,如圖5(a)所示,工況2表示爆炸發(fā)生在艙室3內(nèi)部,如圖5(b)所示。均勻高壓氣云均設(shè)置在艙室中間(見圖5中的淺色區(qū)域),周圍深色區(qū)域表示空氣域,空氣密度為 1.225 kg/m3,壓力為 1.0×105Pa。
本節(jié)主要根據(jù)數(shù)值模擬結(jié)果分析工況1爆炸后爆炸波傳播過程及不同艙室爆炸載荷分布特性。
4.3.1 爆炸波傳播過程分析
為了探討工況1氣云爆炸后爆炸波在多艙室內(nèi)部的傳播過程,這里給出不同時刻爆炸壓力云圖和數(shù)值紋影圖(見圖6)。根據(jù)圖6可知工況1多艙室內(nèi)爆炸過程爆炸波的傳播路徑:爆炸波首先在艙室1中進行自由膨脹運動(圖6(a))爆炸波初次到達艙室壁面時發(fā)生正規(guī)則反射,壁面處形成局部高壓區(qū)域(圖6(b)),之后爆炸波沿著連接導管進入艙室2和艙室3內(nèi)部(圖6(c)),并緊接著在艙室2和艙室3上下壁面處進行系列正規(guī)則斜反射和馬赫反射(圖6(d)和圖6(e)),當爆炸波傳播到艙室2和艙室3的最右端壁面處時再次形成復雜的壁面反射波(圖6(f))。由于多艙室處于封閉空間,因此可以預知爆炸波經(jīng)過上述在不同艙室之間的傳播過程,爆炸波強度將不斷衰弱,爆炸氣體的動能逐漸轉(zhuǎn)化為氣體的內(nèi)能,最終艙室內(nèi)部將趨近于一個準靜態(tài)過程。
4.3.2 爆炸載荷特性分析
為了探討工況1氣云爆炸后不同艙室內(nèi)部爆炸載荷分布特性,這里給出不同艙室內(nèi)部壁面測點超壓時間歷程曲線和沖量歷程時間曲線。針對工況1,由于艙室2和艙室3上下對稱,因此,這里只需要討論艙室1和艙室2內(nèi)部壁面測點的爆炸載荷特性。圖7給出艙室1壁面測點NO1,NO2爆炸載荷時間歷程,從圖7可知艙室內(nèi)爆炸載荷主要呈現(xiàn)峰值不斷衰減的沖擊波載荷和準靜態(tài)壓力。雖然測點NO1,NO2與爆源處于近似相同的距離,但是測點NO1超壓載荷及沖量明顯大于測點NO2,通過分析發(fā)現(xiàn)這主要是由于測點NO2靠近艙室2和艙室3,當爆炸波通過連接導管傳播到艙室2和艙室3內(nèi)部后將明顯減弱測點NO2處的沖擊波強度。
圖5 初始條件及爆炸工況Fig.5 Initial condition and explosion cases
圖6 工況 1 不同時刻爆炸波傳播過程Fig.6 Blast wave propagation process at different time for explosion case one
圖7 艙室 1 壁面測點 NO1,NO2 爆炸載荷時間歷程(工況1)Fig.7 Blast load time histories of gauging point NO1, NO2 in cabin 1 for explosion case one
圖8 給出艙室2壁面測點NO3,NO4爆炸載荷時間歷程,從圖8可知在爆炸初期,測點NO3超壓載荷及沖量稍微大于測點NO4,而在爆炸后期,測點NO3超壓載荷及沖量演變?yōu)樾∮跍y點NO4,通過分析發(fā)現(xiàn)這主要是由于爆炸初期測點NO3相較測點NO4更靠近爆源,因此其接觸到的初始沖擊波強度較大。而在爆炸后期由于測點NO4與連接導管正對,在后續(xù)的爆炸波傳播過程中,其將會受到爆炸波的連續(xù)直接沖擊,反而造成其超壓載荷、沖量時間歷程稍微高于測點NO3。
圖9給出艙室1,艙室2壁面測點NO1,NO2,NO3,NO4爆炸載荷時間歷程,從圖9可知艙室1內(nèi)部測點NO1,NO2爆炸載荷強度明顯強于艙室2內(nèi)部測點NO3,NO4,通過分析發(fā)現(xiàn)這主要是由于爆炸發(fā)生在艙室1內(nèi)部,因此其主要能量集中在艙室1內(nèi)部。同時根據(jù)圖9(b)中的沖量歷程曲線還得知,不同艙室內(nèi)部最終的準靜態(tài)壓力相同,如果將沖量積分曲線開始變?yōu)橹本€段的轉(zhuǎn)折點作為準靜態(tài)壓力開始形成的時間,則在該爆炸工況下不同艙室內(nèi)部測點到達準靜態(tài)壓力的時間基本相同(近似為0.5 ms)。
圖8 艙室 2 壁面測點 NO3,NO4 爆炸載荷時間歷程(工況1)Fig.8 Blast load time histories of gauging point NO3, NO4 in cabin 2 for explosion case one
圖9 艙室 1 及艙室 2 壁面測點 NO1,NO2,NO3,NO4 爆炸載荷時間歷程(工況1)Fig.9 Blast load time histories of gauging point NO1, NO2, NO3,NO4 in cabin 1, 2 for explosion case one
圖10 工況 2 不同時刻爆炸波傳播過程Fig.10 Blast wave propagation process at different time for explosion case two
本節(jié)主要根據(jù)數(shù)值模擬結(jié)果分析工況2爆炸后爆炸波傳播過程及不同艙室爆炸載荷分布特性。
4.4.1 爆炸波傳播過程分析
為了探討工況2氣云爆炸后爆炸波在多艙室內(nèi)部的傳播過程,這里給出不同時刻爆炸壓力云圖和數(shù)值紋影圖(見圖10)。根據(jù)圖10可知工況2多艙室內(nèi)爆炸過程爆炸波的傳播路徑:爆炸波首先在長方形艙室3中進行自由膨脹運動,初次到達艙室上下壁面時發(fā)生正規(guī)則反射,壁面處形成局部高壓區(qū)域(圖10(a)),之后爆炸波沿著連接導管進入艙室1內(nèi)部并以球形方式進行擴張(圖10(b)和圖10(c)),在擴張過程中由于艙室1上壁面的限制,形成局部馬赫反射波(圖10(d)),緊接著,爆炸波傳播到艙室2內(nèi)部并形成復雜的壁面反射波(圖10(e)和圖10(f))。
圖11 艙室 1 壁面測點 NO1,NO2 爆炸載荷時間歷程(工況2)Fig.11 Blast load time histories of gauging point NO1, NO2 in cabin 1 for explosion case two
圖12 艙室 2 壁面測點 NO3,NO4 爆炸載荷時間歷程(工況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)部壁面測點超壓時間歷程曲線和沖量歷程時間曲線。圖11給出艙室1壁面測點NO1,NO2爆炸載荷時間歷程,從圖11可知測點NO1超壓載荷及沖量稍微大于測點NO2,通過分析發(fā)現(xiàn)這主要是由于測點NO1靠近艙室3的連接導管,當爆炸波通過連接導管傳播到艙室1內(nèi)部時,測點NO1將遭受到更強的爆炸波沖擊。
圖12給出艙室2壁面測點NO3,NO4爆炸載荷時間歷程,從圖12可知測點NO3,NO4雖然超壓載荷時間歷程差異較大,但是沖量載荷歷程基本相同,通過分析發(fā)現(xiàn)這主要是由于測點NO3,NO4位于艙室2內(nèi)部,相較于艙室1中的爆源處于遠場區(qū)域,當爆炸波通過連接導管擴張到艙室1再傳播到艙室2內(nèi)部后,由于艙室2容積也較艙室1偏小,傳入的爆炸波沖擊很快使得艙室2內(nèi)部達到準靜態(tài)過程。
圖13給出艙室3壁面測點NO5,NO6爆炸載荷時間歷程,從圖13可知測點NO6初始超壓載荷明顯強于測點NO5,后期強度衰減速度快于測點NO5,這主要是由于艙室3為長方形艙室,測點6距離爆源較近,因此初期沖擊波載荷強度大,由于測點6距離連接導管近,因此其后期超壓衰減速率更快。艙室3容積較小,雖然測點NO6距離連接導管較近,對該處的沖量歷程有一定的影響,但是測點6距離爆源較近,最終使得NO5,NO6沖量歷程曲線基本重合。
圖14給出艙室1~艙室3壁面測點NO1~NO6爆炸載荷時間歷程,從圖14可知艙室3內(nèi)部測點NO5,NO6爆炸載荷強度明顯強于艙室1內(nèi)部測點NO1,NO2,通過分析發(fā)現(xiàn)這主要是由于爆炸發(fā)生在艙室3內(nèi)部,因此其主要能量集中在艙室3內(nèi)部。艙室1內(nèi)部測點NO1,NO2爆炸載荷強度稍微強于艙室2內(nèi)部測點NO3,NO4,通過分析發(fā)現(xiàn)這主要是由于NO1,NO2相較測點NO3,NO4距離爆炸源更近一些,沖擊波強度隨著傳播距離的增大,強度逐漸衰弱。同時根據(jù)圖14(b)中的沖量歷程曲線還得知,盡管爆炸波達到各測點的時間存在差異性,但是該爆炸工況下不同艙室內(nèi)部測點到達準靜態(tài)壓力的時間基本相同(近似為0.5 ms),這一結(jié)論與4.3.2小節(jié)中的分析結(jié)果相同。因此,通過這2個工況的爆炸分析,可以初步得出如下結(jié)論:對于特定的封閉艙室,準靜態(tài)壓力形成時間主要與艙室容積和爆源能量相關(guān)。其實,從準靜態(tài)壓力形成的物理過程來分析也很容易闡述這種關(guān)系,由于準靜態(tài)壓力全場均勻,因此其在全場空間內(nèi)必然應同時到達。
圖13 艙室 3 壁面測點 NO5,NO6 爆炸載荷時間歷程(工況2)Fig.13 Blast load time histories of gauging point NO5, NO6 in cabin 3 for explosion case two
圖14 艙室 1~艙室 3 壁面測點 NO1,NO2,NO3,NO4,NO5,NO6爆炸載荷時間歷程(工況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壁面測點NO1爆炸載荷時間歷程對比圖Fig.15 Comparisons of blast load time histories of gauging point NO1between explosion case one and two
為了進一步分析爆源位置對多艙室內(nèi)部不同測點處爆炸載荷特性的影響規(guī)律,圖15~圖20分別給出2種工況下測點NO1~NO6的超壓載荷和沖量時間歷程曲線圖。
根據(jù)圖15和圖16可知,針對艙室1內(nèi)部測點NO1,NO2,工況1相較工況2具有更強的超壓載荷和更大的沖量,這主要是由于工況1中爆炸發(fā)生在艙室1中。
根據(jù)圖17和圖18可知,針對艙室2內(nèi)部測點NO3,NO4,2種工況超壓載荷歷程差異較大,但是沖量載荷歷程差別較小,這主要是由于艙室2中的測點相較2種工況中的爆源距離較大,屬于遠場的范圍,沖量的歷程主要由準靜態(tài)壓力決定。
根據(jù)圖19和圖20可知,針對艙室3內(nèi)部測點NO5,NO6,工況2相較工況1具有更強的超壓載荷和更大的沖量,這主要是由于工況2中爆炸發(fā)生在艙室3中。
圖16 工況1及工況2壁面測點NO2爆炸載荷時間歷程對比圖Fig.16 Comparisons of blast load time histories of gauging point NO2 between explosion case one and two
圖17 工況1及工況2壁面測點NO3爆炸載荷時間歷程對比圖Fig.17 Comparisons of blast load time histories of gauging point NO3 between explosion case one and two
圖18 工況1及工況2壁面測點NO4爆炸載荷時間歷程對比圖Fig.18 Comparisons of blast load time histories of gauging point NO4 between explosion case one and two
圖19 工況1及工況2壁面測點NO5爆炸載荷時間歷程對比圖Fig.19 Comparisons of blast load time histories of gauging point NO5 between explosion case one and two
圖20 工況1及工況2壁面測點NO5爆炸載荷時間歷程對比圖Fig.20 Comparisons of blast load time histories of gauging point NO5 between explosion case one and two
采用三階WENO有限差分格式,自主開發(fā)了密閉多艙室內(nèi)爆炸波二維數(shù)值計算程序。利用所開發(fā)的程序開展了多艙室內(nèi)柱形氣云爆炸波傳播過程及不同艙室內(nèi)部爆炸載荷特性研究。通過本文的研究主要得到以下幾點結(jié)論:
1)多艙室內(nèi)爆炸過程中,爆炸所在艙室內(nèi)部爆炸載荷強度明顯強于相鄰艙室,尤其需要重點關(guān)注爆炸所在艙室通過連接導管對直接相鄰艙室爆炸載荷的影響;爆源位置對不同艙室內(nèi)部超壓載荷時間歷程影響較大,而對處于遠場的艙室內(nèi)沖量時間歷程影響較小。
2)對于特定的封閉艙室,準靜態(tài)壓力峰值及準靜態(tài)壓力形成時間均全場均勻,準靜態(tài)壓力形成時間主要與艙室容積和爆源能量相關(guān)。
3)本文研究可為工程抗爆結(jié)構(gòu)設(shè)計及爆炸毀傷評估提供一定的參考和指導。