固體火箭發(fā)動(dòng)機(jī)結(jié)構(gòu)意外爆炸的數(shù)值模擬與風(fēng)險(xiǎn)評估?

樊凱旋 彭 旭 任家帆 饒國寧

南京理工大學(xué)化工學(xué)院(江蘇南京，210094)

引言

近年來，隨著我國航天航空技術(shù)的不斷快速發(fā)展，對于航天航空產(chǎn)品作業(yè)過程中的安全要求越來越高。 作為航天航空產(chǎn)品的動(dòng)力裝置，固體火箭發(fā)動(dòng)機(jī)是不利用周圍介質(zhì)，而是利用自身攜帶的固體推進(jìn)劑來生成工質(zhì)的直接反作用發(fā)動(dòng)機(jī)。 固體推進(jìn)劑是化學(xué)推進(jìn)劑的一種，它利用化學(xué)反應(yīng)所釋放出的熱能作為推進(jìn)系統(tǒng)的熱源[1]。 隨著現(xiàn)代固體推進(jìn)劑能量的提高和火炸藥技術(shù)的融合，固體推進(jìn)劑中會(huì)含有高能炸藥成分，在外界刺激作用下一旦發(fā)生爆炸，將會(huì)導(dǎo)致災(zāi)難性事故的發(fā)生[2-4]。

從20 世紀(jì)60 年代起，國外對固體火箭發(fā)動(dòng)機(jī)的意外爆炸及其后果開展了深入的研究[5-7]，不僅進(jìn)行了大量的模擬仿真和理論分析工作，還進(jìn)行過多次的固體發(fā)動(dòng)機(jī)爆炸試驗(yàn)，積累了大量的試驗(yàn)數(shù)據(jù)，并且依此制定了安全防范措施。目前，國內(nèi)的研究工作主要集中在幾個(gè)方面：1)針對固體發(fā)動(dòng)機(jī)在外來沖擊下是否安全[8]，并未涉及到爆炸所產(chǎn)生的危害；2)固體火箭發(fā)動(dòng)機(jī)自毀的危險(xiǎn)性研究[9-10]；3)高能固體火箭發(fā)動(dòng)機(jī)爆炸沖擊波毀傷效應(yīng)研究[11]，這里只是涉及在無遮蔽的情況下固體發(fā)動(dòng)機(jī)爆炸沖擊波對人的殺傷作用，并沒有考慮到固體發(fā)動(dòng)機(jī)在工房中爆炸的情況。 而針對航天企業(yè)裝配作業(yè)時(shí)工房內(nèi)的固體發(fā)動(dòng)機(jī)爆炸事故的模擬和后果評估方面鮮有研究報(bào)道。

本文中，以航天產(chǎn)品裝配工房作業(yè)為研究對象，采用顯式動(dòng)力學(xué)軟件AUTODYN 模擬固體發(fā)動(dòng)機(jī)作業(yè)過程中的意外爆炸，分析工房內(nèi)的沖擊波流場分布規(guī)律以及對人員傷害的影響，研究結(jié)果可以為裝配作業(yè)的風(fēng)險(xiǎn)評估、安全對策與防范提供理論支撐。

1 有限元模型的建立

以某型號(hào)固體火箭發(fā)動(dòng)機(jī)裝配作業(yè)為研究對象，只考慮沖擊波對人的傷害作用，所以可以將固體發(fā)動(dòng)機(jī)等效為相應(yīng)質(zhì)量的TNT。 有限元模型如圖1所示。

1.1 材料方程及參數(shù)的確定

1)空氣。 空氣采用Ideal Gas 理想氣體狀態(tài)方程，其表達(dá)式為

式中：p為壓力；ρ為空氣密度，ρ ＝1.225 kg/m3；e為氣體單位質(zhì)量內(nèi)能，e ＝206.8 kJ/g；γ為絕熱指數(shù)，γ ＝1.4。

2)炸藥。 根據(jù)固體發(fā)動(dòng)機(jī)的TNT 當(dāng)量，采用JWL 狀態(tài)方程來描述其爆轟過程，具體表達(dá)式為

式中：p為爆壓；V為爆轟產(chǎn)物的相對體積；E為初始比內(nèi)能；A1、B1、R1、R2和ω為材料常數(shù)，其具體參數(shù)取值如表1 所示。 表1 中，ρ為TNT 的密度；D為炸藥的爆速[12-13]。

3)水泥混凝土。 CONC-35MPA 混凝土的狀態(tài)方程為p-alpha，其表達(dá)式為

表1 炸藥的JWL 狀態(tài)方程參數(shù)Tab.1 Parameters of JWL equation of state for explosives

式中：ρ0為材料初始密度；ρ為壓縮過程中材料的密度；e為初始內(nèi)能；A1、A2、A3、B0、B1、T1、T2為材料常數(shù)，其具體參數(shù)取值如表2 所示[14]。

表2 CONC-35MPA 混凝土的材料方程參數(shù)Tab.2 Material equation parameters of CONC-35MPA concrete

強(qiáng)度模型則選擇適用于描述爆炸過程大變形、高應(yīng)變的RHT Concrete 本構(gòu)模型。 RHT 本構(gòu)模型引入了3 個(gè)極限面， 即彈性極限面、失效面和殘余強(qiáng)度面，分別描述混凝土的初始屈服強(qiáng)度、失效強(qiáng)度及殘余強(qiáng)度的變化規(guī)律。 失效面的方程為

表3 RHT 本構(gòu)模型強(qiáng)度參數(shù)Tab.3 Strength parameters of RHT constitutive model

1.2 數(shù)值計(jì)算模型

AUTODYN 里建立的模型為空氣模型、產(chǎn)品模型與工房模型。 氣域?yàn)殚L方體結(jié)構(gòu)，將產(chǎn)品和工房包括在內(nèi)，將三者用完全耦合的方式連接，空氣域的邊界條件為Flow-out (無限域、無反射邊界條件)。

1.2.1 空氣模型

空氣模型采用Euler，3D-Multi-Material 算法，空氣域?yàn)?00 m×100 m×15 m 的長方體結(jié)構(gòu)，采用六面體結(jié)構(gòu)化網(wǎng)格，將網(wǎng)格設(shè)置為漸變網(wǎng)格，最小網(wǎng)格尺寸為50 mm。 網(wǎng)格劃分如圖2 所示。

1.2.2 產(chǎn)品模型

根據(jù)《軍工燃燒爆炸品工程安全技術(shù)規(guī)范征求意見稿2019》可以得出，此型號(hào)固體火箭發(fā)動(dòng)機(jī)的TNT 當(dāng)量為0.4[17]，該數(shù)值是通過實(shí)驗(yàn)得出的。

此產(chǎn)品的TNT 當(dāng)量比為0. 4，TNT 的密度為1.63 g/cm3，密度體積公式為

根據(jù)式(7)可以得到TNT 的尺寸為0.55 m ×0.40 m ×0.40 m。 TNT 采用Euler 算法，單發(fā)產(chǎn)品炸點(diǎn)為炸藥的中心，雙發(fā)產(chǎn)品的炸點(diǎn)有兩個(gè)，分別布置在兩個(gè)炸藥的中心，起爆時(shí)間為同時(shí)起爆，雙發(fā)產(chǎn)品間距為1 m，炸藥嵌在空氣域里。 炸藥采用的是六面體結(jié)構(gòu)化網(wǎng)格，網(wǎng)格尺寸為50 mm，單發(fā)和雙發(fā)產(chǎn)品網(wǎng)格劃分如圖3 和圖4 所示。

1.2.3 工房模型

工房為長方體狀，尺寸為30 m×15 m×9 m，工房采用六面體結(jié)構(gòu)化網(wǎng)格，網(wǎng)格尺寸為200 mm。 單發(fā)產(chǎn)品放置在工房地面的中心處，雙發(fā)產(chǎn)品在工房地面的中軸線上并排放置。 網(wǎng)格劃分如圖5 所示。

為了定性和定量地分析工房內(nèi)部壓力分布規(guī)律，在空氣域內(nèi)設(shè)置了一系列監(jiān)測點(diǎn)，監(jiān)測點(diǎn)的位置均為距工房地面1 m 高。 為了便于比較，設(shè)置單發(fā)產(chǎn)品和雙發(fā)產(chǎn)品爆炸時(shí)監(jiān)測點(diǎn)的位置相同，如圖6所示。

2 數(shù)值模擬結(jié)果分析

產(chǎn)品在工房地面爆炸時(shí)，形成的爆炸沖擊波向四周傳播，當(dāng)遇到壁面時(shí)反射形成反射沖擊波，這些反射波及入射波相互疊加，使空間內(nèi)的壓力增大。TNT 在密閉空間爆轟后，能量的釋放主要分為兩個(gè)過程：爆轟和后燃燒過程。 爆轟形成了沖擊波壓力，后燃燒形成了密閉空間中的準(zhǔn)靜態(tài)壓力，高頻沖擊波壓力波與低頻準(zhǔn)靜態(tài)壓力波疊加[18]。 這里只研究爆轟產(chǎn)生的沖擊波超壓的影響，以此來進(jìn)行定性和定量分析。

2.1 不同時(shí)刻的壓力分布

單發(fā)產(chǎn)品爆炸過程中，不同時(shí)刻工房的壓力云如圖7 所示。

當(dāng)t＝0 ms 時(shí)，由于爆炸還沒有開始，工房處在空氣域內(nèi)，所以此時(shí)工房內(nèi)的壓力應(yīng)該與環(huán)境壓力相同，壓力為101 kPa。

當(dāng)t為10 ms時(shí)，窗口處的壓力大于工房內(nèi)其他位置的壓力，窗口處的最大壓力為1 018 kPa。這是由于：炸藥爆炸沖擊波以球形的方式向外傳播，炸藥距離窗口比距離墻面更近，所以沖擊波到達(dá)得更早，窗口的壓力要比其余位置大；泄壓時(shí)，工房內(nèi)部的壓力從窗口泄放，就會(huì)導(dǎo)致窗口附近的壓力急劇升高。

當(dāng)t為20 ms 時(shí)，工房頂部的最大壓力可以達(dá)到2 800 kPa，工房內(nèi)大部分區(qū)域的壓力為負(fù)值。

當(dāng)t為30 ms 時(shí)，工房內(nèi)的壓力值維持在450 kPa 左右，工房頂部的壓力開始下降。

當(dāng)t為40 ms 時(shí)，工房內(nèi)大部分區(qū)域壓力在620 kPa 左右，最大壓力為5 600 kPa。

t從50 ms 到70 ms，窗口處的壓力為工房內(nèi)的最大壓力，最大壓力從5 700 kPa 上升到8 000 kPa，工房內(nèi)出現(xiàn)負(fù)壓區(qū)域。

雙發(fā)產(chǎn)品爆炸過程中，不同時(shí)刻工房的壓力云如圖8 所示。

模擬未開始時(shí)，整個(gè)工房的壓力與環(huán)境壓力相同；當(dāng)t為10 ms 時(shí)，工房頂部和窗口處的壓力率先到達(dá)最大值，最大壓力為2 600 kPa；當(dāng)t為20 ms時(shí)，工房頂部的最大壓力可以達(dá)到3 380 kPa，工房內(nèi)大部分區(qū)域的壓力為負(fù)值；當(dāng)t為30 ms 時(shí)，工房內(nèi)的壓力為580 kPa，工房頂部的壓力開始下降；當(dāng)t為40 ms 時(shí)，窗口處的壓力依然是工房內(nèi)的最大壓力，最大壓力為9 100 kPa；t從50 ms 到70 ms，窗口處的壓力為工房內(nèi)的最大壓力，最大壓力從10 730 kPa 上升到13 200 kPa。

2.2 沖擊波超壓與傷害半徑

在整個(gè)空氣域內(nèi)要對沖擊波超壓進(jìn)行定量分析，就必須在空氣域的不同位置添加監(jiān)測點(diǎn)，得到實(shí)時(shí)壓力變化情況，根據(jù)壓力變化曲線來判斷沖擊波超壓的大小，再根據(jù)超壓對人的傷害作用，從而求出傷害半徑。 考慮窗口對壓力的變化會(huì)有影響，所以分兩種情況來對單發(fā)產(chǎn)品和雙發(fā)產(chǎn)品爆炸的傷害范圍進(jìn)行分析。 其中，X方向與窗口平行，Y方向與窗垂直。 查閱文獻(xiàn)得到?jīng)_擊波超壓對人的傷害作用如表4 所示[19]。

2.2.1X方向上的壓力分析

表4 沖擊波超壓對人的傷害作用Tab.4 Damage effect of shock wave overpressure on humans

X方向上監(jiān)測點(diǎn)的分布為：沿著工房的中心向X軸分布，所有的點(diǎn)均在與X軸平行的直線上。 測點(diǎn)的位置如表5 所示。

表5 X 方向上測點(diǎn)的位置Tab.5 Location of points on the X direction

各測點(diǎn)的壓力曲線見圖9、圖10。

可以看出，單發(fā)產(chǎn)品爆炸時(shí)，測點(diǎn)1?！?＃的壓力均從101 kPa開始變化，不同位置的沖擊波超壓到達(dá)時(shí)間也不相同，到炸點(diǎn)的距離越遠(yuǎn)，超壓到達(dá)的時(shí)間也越大。 并且隨著到炸點(diǎn)距離的增大，沖擊波超壓也越來越小。 這也符合沖擊波超壓作用的原則。

雙發(fā)產(chǎn)品爆炸時(shí)，測點(diǎn)1＃～5＃的壓力也是從環(huán)境壓力開始變化的，壓力變化趨勢和單發(fā)產(chǎn)品爆炸時(shí)的壓力變化基本相同，但是超壓峰值的增長與爆炸產(chǎn)品數(shù)量的增長呈現(xiàn)非線性關(guān)系。 根據(jù)表4 沖擊波超壓對人的傷害作用，得到單發(fā)產(chǎn)品和雙發(fā)產(chǎn)品爆炸的X方向上的傷害范圍，如表6 所示。

表6 X 方向上的傷害半徑Tab.6 Damage radius in the X direction m

所以，單發(fā)產(chǎn)品爆炸在X方向上無傷害半徑為40.0 m 之外；輕傷半徑＞29.0 ～40.0 m；重傷半徑20.5 ～29.0 m；死亡半徑為20.5 m 之內(nèi)。 雙發(fā)產(chǎn)品爆炸在X方向上無傷害半徑為45.0 m 之外；輕傷半徑＞32.0 ～45.0 m；重傷半徑24.0 ～32.0 m；死亡半徑為24.0 m 之內(nèi)。

2.2.2Y方向上的壓力分析

Y方向上監(jiān)測點(diǎn)的分布為：沿著工房的中心向Y軸分布，所有測點(diǎn)均在與Y軸平行的直線上。 測點(diǎn)的位置如表7 所示。 各測點(diǎn)的壓力曲線如圖11、圖12 所示。

表7 Y 方向上點(diǎn)的位置Tab.7 Location of the point in the Y direction

可以看出，單發(fā)產(chǎn)品爆炸時(shí)，測點(diǎn)7?！?2＃的壓力也是從環(huán)境壓力開始變化。 為了更好地比較各個(gè)點(diǎn)的沖擊波超壓的差異以及壓力變化規(guī)律，這里將測點(diǎn)6＃的壓力曲線舍去，因?yàn)闇y點(diǎn)6＃的沖擊波超壓過大，如果放在同一個(gè)坐標(biāo)系中，很難根據(jù)壓力曲線對其他點(diǎn)做出比較。

圖11 反映出隨著到炸點(diǎn)的距離變遠(yuǎn)，沖擊波超壓到達(dá)的時(shí)間隨之增大，沖擊波超壓也在逐漸減小，同樣符合沖擊波超壓作用的原則。

雙發(fā)產(chǎn)品爆炸時(shí)，測點(diǎn)7?！?2＃壓力增加的幅度明顯要比測點(diǎn)1?！?＃的壓力增加幅度大，這是因?yàn)闇y點(diǎn)7?！?2＃分布于兩個(gè)炸藥連線的中垂線上，由于沖擊波是以球形的方式向外擴(kuò)散，所以測點(diǎn)7?！?2＃受到兩股沖擊波同時(shí)作用。 根據(jù)表4 沖擊波超壓對人的傷害作用，就可以得出在Y方向上的傷害半徑，如表8 所示。

表8 Y 方向上的傷害半徑Tab.8 Damage radius in the Y direction m

所以單發(fā)產(chǎn)品爆炸在Y方向上無傷害半徑為43.0 m 之外，輕傷半徑＞31.0 ～43.0 m；重傷半徑22.0 ～31.0 m；死亡半徑為22.0 m 之內(nèi)。 雙發(fā)產(chǎn)品爆炸在Y方向上無傷害半徑為55.0 m 之外，輕傷半徑＞42.5 ～55.0 m；重傷半徑32.0 ～42.5 m；死亡半徑在32.0 m 之內(nèi)。

2.3 模擬結(jié)果與實(shí)驗(yàn)結(jié)果的比較

根據(jù)陽建紅等[20]做的無遮蔽情況下高能推進(jìn)劑爆炸實(shí)驗(yàn)，可以將實(shí)驗(yàn)結(jié)果與模擬結(jié)果作比較，從而印證上述模擬方法是正確的。

具體做法為：由實(shí)驗(yàn)中推進(jìn)劑的質(zhì)量和推進(jìn)劑TNT 當(dāng)量可以計(jì)算出TNT 質(zhì)量。 由于雙發(fā)產(chǎn)品之間的沖擊波會(huì)相互影響，這里只對單發(fā)產(chǎn)品進(jìn)行分析。 將文章中有限元模型的工房部分刪除，產(chǎn)品質(zhì)量與實(shí)驗(yàn)中換算得出的TNT 質(zhì)量相同，然后再增加測點(diǎn)，模擬測點(diǎn)位置與試驗(yàn)測點(diǎn)位置相同，其余條件保持不變。 可以得出測點(diǎn)處沖擊波超壓的數(shù)據(jù)，如表9 所示。

從表9 可以看出，在距離爆源大于6 m 時(shí)，數(shù)值模擬結(jié)果與實(shí)驗(yàn)所測數(shù)據(jù)誤差范圍在10%以內(nèi)；距離爆源越近，爆炸場內(nèi)沖擊波傳播越復(fù)雜，馬赫反射的影響越顯著，從而造成了模擬數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)相差較大。 通過實(shí)驗(yàn)結(jié)果與模擬結(jié)果相比較，說明文中所建立的計(jì)算模型與方法是正確的。

3 結(jié)論

通過數(shù)值模擬的手段研究了航天產(chǎn)品裝配工房內(nèi)單發(fā)產(chǎn)品和雙發(fā)產(chǎn)品的爆炸過程，得出如下結(jié)論：

表9 實(shí)驗(yàn)結(jié)果與模擬結(jié)果相比較Tab.9 Result comparison of experiment with simulation

1)根據(jù)不同時(shí)刻的壓力云圖可以得到，產(chǎn)品在工房內(nèi)爆炸要比在空氣中爆炸復(fù)雜得多，超壓峰值的大小與工房的結(jié)構(gòu)有很大的關(guān)系。

2)當(dāng)單發(fā)產(chǎn)品爆炸時(shí)，沖擊波超壓在X方向上對人的殺傷距離為29.0 m，爆炸后人員的安全距離為40.0 m；在Y方向上對人的殺傷距離為31.0 m，爆炸后人員的安全距離為43.0 m。 當(dāng)雙發(fā)產(chǎn)品爆炸時(shí)，沖擊波超壓在X方向上對人的殺傷距離為32.0 m，爆炸后人員的安全距離為45.0 m；在Y方向上對人的殺傷距離為42.5 m，爆炸后人員的安全距離為55.0 m。