航天發(fā)射場液體推進(jìn)劑的泄漏擴(kuò)散模型研究*

柳寧遠(yuǎn)，崔村燕，辛騰達(dá)，欒　驍

(1. 航天工程大學(xué) 研究生管理大隊(duì)，北京 101416；2. 航天工程大學(xué) 航天裝備系，北京 101416；3. 中國人民解放軍63618部隊(duì)，新疆 庫爾勒 841001 )

0　引言

隨著航天事業(yè)的發(fā)展，發(fā)射任務(wù)密集，發(fā)射場安全問題顯得尤為重要。液體導(dǎo)彈推進(jìn)劑UDMH毒性較強(qiáng)，一旦發(fā)生爆炸或推進(jìn)劑泄漏事故，有毒推進(jìn)劑的逸散會(huì)對(duì)發(fā)射塔架周邊的人員安全和環(huán)境造成巨大危害[1]。因此，有必要對(duì)泄漏事故進(jìn)行快速評(píng)估，分析不同泄漏時(shí)間有毒推進(jìn)劑在大氣中的擴(kuò)散范圍以及人員危害區(qū)域，以便及時(shí)采取防護(hù)措施，保證人員的安全。

對(duì)危害性氣體擴(kuò)散過程的數(shù)學(xué)模擬以及正確估算各種條件下危害性氣體濃度隨時(shí)空的分布和變化，關(guān)鍵在于擴(kuò)散模型在各種條件下的應(yīng)用。國內(nèi)外關(guān)于不同條件下氣體擴(kuò)散模型的應(yīng)用有比較多的研究，包括高架源氣體泄漏，濱海地形大氣擴(kuò)散，城市污染氣體擴(kuò)散，火災(zāi)煙霧擴(kuò)散等[2];關(guān)于發(fā)射場推進(jìn)劑泄漏情形下的擴(kuò)散研究較少，陳新華等提出了爆炸后殘余推進(jìn)劑蒸發(fā)形成的毒氣在大氣中擴(kuò)散的解析解[3]，但對(duì)于擴(kuò)散模型在該情形下的適用性有待進(jìn)一步研究。

本文針對(duì)擴(kuò)散模型在航天發(fā)射場推進(jìn)劑氣體擴(kuò)散的適用性進(jìn)行研究，結(jié)合推進(jìn)劑偏二甲肼在航天發(fā)射場蒸發(fā)擴(kuò)散的條件，對(duì)擴(kuò)散模型進(jìn)行完善并建立數(shù)學(xué)模型，并通過擴(kuò)散模型計(jì)算結(jié)果和對(duì)應(yīng)數(shù)值模擬結(jié)果及已有實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比分析，為航天發(fā)射場擴(kuò)散模型的建立和完善提供參考。擴(kuò)散模型的建立有利于事故后快速估算有毒推進(jìn)劑濃度分布，為泄漏事故危害性評(píng)估和防護(hù)提供依據(jù)。

1　擴(kuò)散模型建立

1.1　理論模型及基本假設(shè)

根據(jù)梯度輸送理論，由湍流運(yùn)動(dòng)引起的局部位置質(zhì)量通量與該位置擴(kuò)散物質(zhì)的平均濃度梯度成正比。式(1)為根據(jù)梯度輸送理論導(dǎo)出的普遍形式湍流擴(kuò)散方程，說明流體中某物質(zhì)的散布是由湍流擴(kuò)散所引起的。x，y，z方向的湍流擴(kuò)散系數(shù)分別為Kx，Ky，Kz，擴(kuò)散系數(shù)為時(shí)空的函數(shù)，與流場形式有關(guān)[4]。

(1)

若將坐標(biāo)系x軸與風(fēng)向取一致，z軸垂直向上，假設(shè)流場在3個(gè)方向的擴(kuò)散系數(shù)為常數(shù)(即斐克擴(kuò)散情形)，則湍流擴(kuò)散方程可簡化為：

(2)

(3)

為簡化計(jì)算，略去x方向的湍流項(xiàng)。

當(dāng)x>0時(shí)，

(4)

當(dāng)z→0時(shí)，

(5)

求解可得，有風(fēng)時(shí)連續(xù)點(diǎn)源濃度分布為：

(6)

(7)

1.2　擴(kuò)散模型在發(fā)射場推進(jìn)劑蒸發(fā)擴(kuò)散中的應(yīng)用

航天發(fā)射場附近地勢平坦開闊，氣象穩(wěn)定，符合平穩(wěn)均勻湍流的假設(shè)，粒子擴(kuò)散位移符合正態(tài)分布形式，風(fēng)速一般大于1 m/s，符合上節(jié)所述擴(kuò)散模型和基本假設(shè)，但根據(jù)實(shí)際情況需要進(jìn)行改進(jìn)完善。

液體火箭意外發(fā)生爆炸時(shí)，部分推進(jìn)劑在高溫下蒸發(fā)成高溫氣體迅速擴(kuò)散，部分液體推進(jìn)劑散布在周邊緩慢蒸發(fā)。計(jì)算模型只適用于緩慢蒸發(fā)的推進(jìn)劑，高溫氣體隨沖擊波迅速擴(kuò)散至遠(yuǎn)處，相對(duì)影響較小，在計(jì)算中不予以考慮。

1.2.1殘留推進(jìn)劑泄漏源強(qiáng)

殘留推進(jìn)劑散布在不規(guī)則區(qū)域，面積為Aw，對(duì)于面源的源強(qiáng)計(jì)算，可通過污染區(qū)面積和污染區(qū)推進(jìn)劑蒸發(fā)速率來估計(jì)：

(9)

以偏二甲肼(UDMH)為例，偏二甲肼在土壤中的蒸發(fā)速率為[5]：

(10)

式中：v為大氣環(huán)境下土壤中液體蒸發(fā)速率，kg/(m2·s)；ps,f為UDMH飽和蒸汽壓，Pa；p為發(fā)射場當(dāng)?shù)卮髿鈮海琍a；u∞為平均風(fēng)速，m/s；up為參考風(fēng)速，計(jì)算時(shí)取5m/s。

液體推進(jìn)劑偏二甲肼與四氧化二氮分子量比空氣大，且易于反應(yīng)分解。對(duì)于氣體和很小的粒子沉降作用可忽略不考慮，但是由于湍流擴(kuò)散和布朗運(yùn)動(dòng)沉積到各種表面，被吸收或者發(fā)生其他反應(yīng)會(huì)對(duì)其有一定消耗，故為對(duì)公式中的泄漏源強(qiáng)度Q進(jìn)行修正，以考慮空氣中污染物質(zhì)的量的減損[6]。

(11)

式中：vd為沉積速度，m/s，取值0.01。

1.2.2地面反射效應(yīng)

殘留推進(jìn)劑分布在地表，屬于地面泄漏源。由于地面的反射造成污染物濃度的上升，為簡化計(jì)算，假設(shè)為全反射的情況，泄漏濃度則是無界條件下的2倍，即：

(12)

1.2.3殘留推進(jìn)劑等效面源擴(kuò)散

呈面塊散布的污染物排放源稱為面源，按照點(diǎn)源擴(kuò)散公式可以沿x，y方向積分得到。

(13)

為了簡化模型快速計(jì)算，也可以采用虛點(diǎn)源法，將面源化為點(diǎn)源處理。使得由虛點(diǎn)源排放的污染物經(jīng)過虛擬距離x0后，與面源具有同樣的擴(kuò)散幅。

(14)

常采用經(jīng)驗(yàn)方法給出初始擴(kuò)散幅：

(15)

則地面濃度有：

(16)

空間濃度分布：

(17)

1.2.4擴(kuò)散參數(shù)的確定方法

根據(jù)湍流統(tǒng)計(jì)理論，氣體在隨機(jī)流場中的擴(kuò)散能力和散布范圍由大氣擴(kuò)散參數(shù)來體現(xiàn)。對(duì)應(yīng)x，y，z3個(gè)方向分別用σx，σy，σz表示，隨著下風(fēng)距離x變化而變化。當(dāng)風(fēng)速大于一定值時(shí)，x方向的湍流擴(kuò)散可以忽略不計(jì)。所以，在不同的高斯煙羽模式中，確定大氣擴(kuò)散參數(shù)的關(guān)鍵在于確定不同條件下不同的擴(kuò)散參數(shù)σy，σz。擴(kuò)散系數(shù)σy，σz的大小與大氣湍流結(jié)構(gòu)，離地高度，地面粗糙度，泄漏持續(xù)時(shí)間，抽樣時(shí)間間隔，風(fēng)速，以及離泄漏源的距離等因素有關(guān)。根據(jù)地面氣象觀測數(shù)據(jù)，對(duì)大氣的擴(kuò)散能力進(jìn)行判別，分為A，B，C，D，E，F(xiàn)共6類大氣穩(wěn)定度擴(kuò)散級(jí)別。每個(gè)擴(kuò)散參數(shù)對(duì)應(yīng)每1類大氣穩(wěn)定度有1條擴(kuò)散曲線，表示擴(kuò)散參數(shù)隨下風(fēng)距離x的變化。根據(jù)曲線就可以確定某區(qū)域確定穩(wěn)定度下下風(fēng)距離x的大氣擴(kuò)散參數(shù)的值。

1)確定大氣穩(wěn)定度級(jí)別

大氣穩(wěn)定度可以分為A，B，C，D，E，F(xiàn) 6類，按照從不穩(wěn)定到穩(wěn)定依次分級(jí)，如表1～2所示。這里采用的是D.B. Turner所提出的用太陽高度角來定量判定日照強(qiáng)度的方法[7]。

表1　由太陽高度角確定日照強(qiáng)度與等級(jí)

表2　Turner的穩(wěn)定度分級(jí)方法

Turner的分級(jí)方法定量相對(duì)比較確切，只要有地面風(fēng)速，云量和云高的觀測數(shù)據(jù)，就可以確定穩(wěn)定度級(jí)別，是實(shí)際應(yīng)用中更為普遍使用的方法。

2)根據(jù)大氣穩(wěn)定度級(jí)別選擇擴(kuò)散參數(shù)曲線

按照相應(yīng)的大氣穩(wěn)定度級(jí)別確定擴(kuò)散參數(shù)曲線，大氣擴(kuò)散參數(shù)曲線可以有多種表達(dá)方式。比如，參數(shù)表達(dá)式，觀測點(diǎn)數(shù)據(jù)，曲線表示等，可以根據(jù)不同需要按照下風(fēng)向距離x確定擴(kuò)散參數(shù)確定具體位置的大氣擴(kuò)散參數(shù)的值[7]，如表3所示。

表3　不同大氣穩(wěn)定度下擴(kuò)散參數(shù)關(guān)于下風(fēng)向距離x的表達(dá)式

2　數(shù)值模擬

結(jié)合航天發(fā)射場的氣體泄漏的實(shí)際情況，通過理論分析得到的氣體計(jì)算方程，其準(zhǔn)確性需要試驗(yàn)數(shù)據(jù)的驗(yàn)證。目前，比較可行的試驗(yàn)方法有數(shù)值模擬方法和實(shí)測試驗(yàn)方法，數(shù)值模擬計(jì)算全面精確，便于分析，但為簡化條件和實(shí)際情形有一定差異；實(shí)測試驗(yàn)結(jié)果相對(duì)有說服力，但存在測量誤差和偶然誤差。綜合考慮，擬采用對(duì)比分析的方法彌補(bǔ)各自局限性。

2.1　數(shù)值模擬

對(duì)地面源泄漏過程進(jìn)行數(shù)值模擬，假設(shè)泄漏源為連續(xù)泄漏源，泄漏源等效面積為1 m2，位于計(jì)算域下風(fēng)向10 m處地面。計(jì)算區(qū)域高度在100 m以下，不考慮大氣的分層現(xiàn)象，假設(shè)初始條件計(jì)算區(qū)域空氣均勻分布，計(jì)算域的幾何模型如圖1所示。

圖1　幾何模型示意Fig.1　Geometric model

偏二甲肼蒸發(fā)速率主要受到蒸發(fā)面積、風(fēng)速、溫度和濕度等因素的影響，且滿足如下計(jì)算公式：

(18)

式中：Aw為大氣環(huán)境下土壤中液體蒸發(fā)面積，m2；ps,f為UDMH飽和蒸汽壓，Pa，取20 kPa；p為發(fā)射場當(dāng)?shù)卮髿鈮?，Pa；u∞為平均風(fēng)速，m/s；up為參考風(fēng)速，計(jì)算時(shí)取5 m/s。假設(shè)計(jì)算區(qū)域風(fēng)向不變，風(fēng)速為3 m/s，壓強(qiáng)為89 200 Pa，空氣濕度20%，溫度30 ℃。則蒸發(fā)速率為：

(19)

將蒸發(fā)速率作為液池表面的速度邊界條件代入數(shù)值計(jì)算，模擬地面殘留推進(jìn)劑蒸發(fā)的擴(kuò)散。

在數(shù)值計(jì)算中主要采用有限元方法，描述偏二甲肼在開放空間內(nèi)泄漏擴(kuò)散過程的基本控制方程主要有氣體狀態(tài)方程、流體連續(xù)性方程，質(zhì)量方程、動(dòng)量方程、能量方程和組分輸運(yùn)方程。

經(jīng)數(shù)值模擬計(jì)算所得的地面偏二甲肼濃度分布如圖2～3所示。從圖中可以定性分析得到：擴(kuò)散濃度呈對(duì)稱分布，說明擴(kuò)散作用垂直于下風(fēng)向的方向上大小相同。

圖2　600s地面偏二甲肼體積濃度分布Fig.2　600s UDMH volume fraction on the ground

圖3　垂直下風(fēng)向截面偏二甲肼體積濃度分布Fig.3　UDMH volume fraction of vertical downward direction

圖4　下風(fēng)向不同距離地面體積濃度分布Fig.4　Volume fraction curve of different downwind distance on the ground

如圖2～3所示，擴(kuò)散的偏二甲肼氣體主要分布在下風(fēng)向，下風(fēng)向分布范圍遠(yuǎn)大于其他方向，說明平流輸送對(duì)擴(kuò)散的作用遠(yuǎn)大于湍流擴(kuò)散的作用。

取下風(fēng)向10，30，50，90 m處，地面偏二甲肼濃度繪制曲線，如圖4所示。濃度分布符合正態(tài)分布，與計(jì)算公式及統(tǒng)計(jì)規(guī)律相吻合。

2.2　實(shí)驗(yàn)結(jié)果對(duì)比驗(yàn)證

2.2.1實(shí)驗(yàn)條件

推進(jìn)劑擴(kuò)散實(shí)驗(yàn)危險(xiǎn)性較大，成本較高[11]，故采用已有實(shí)驗(yàn)數(shù)據(jù)進(jìn)行相關(guān)比較與驗(yàn)證。陳新華等針對(duì)液體火箭推進(jìn)劑爆炸毒氣逸散做過大量的實(shí)驗(yàn)研究[8]，以其中相關(guān)實(shí)驗(yàn)結(jié)果與數(shù)值模擬及計(jì)算公式結(jié)果進(jìn)行對(duì)比。

實(shí)驗(yàn)?zāi)M推進(jìn)劑爆炸后毒氣擴(kuò)散過程，大氣中毒氣濃度分布情況，爆炸后計(jì)算地面殘留推進(jìn)劑量，并對(duì)地面殘留推進(jìn)劑蒸發(fā)擴(kuò)散的濃度分布進(jìn)行檢測分析。大氣中濃度采樣是利用電動(dòng)吸氣式大氣采樣器進(jìn)行，9個(gè)采樣點(diǎn)分布在下風(fēng)向不同距離，有毒氣體通過吸氣管進(jìn)入吸收液，事后對(duì)吸收液進(jìn)行檢測分析。

實(shí)驗(yàn)的主要環(huán)境條件[9]為：氣壓89 200 Pa，溫度29℃，風(fēng)速1～3 m/s，空氣濕度19%，大氣穩(wěn)定度D類。環(huán)境條件與數(shù)值模擬的假設(shè)條件基本一致，可以認(rèn)為是相同條件下的蒸發(fā)擴(kuò)散過程。

2.2.2采樣點(diǎn)數(shù)據(jù)比較

采樣點(diǎn)位置的數(shù)值模擬結(jié)果和擴(kuò)散模型計(jì)算結(jié)果如表4所示。

表4　實(shí)驗(yàn)測量值、數(shù)值模擬值和擴(kuò)散模型值

圖5　采樣點(diǎn)位置實(shí)驗(yàn)值Q1、數(shù)值模擬值Q2和擴(kuò)散模型計(jì)算值Q3比較Fig.5　Comparison between test measurements,numerical simulations and diffusion model values

比較數(shù)值模擬值與擴(kuò)散模型計(jì)算值，擴(kuò)散模型計(jì)算值整體小于數(shù)值模擬值，兩者相關(guān)度較大，同一位置濃度偏差在15%左右，最大不超過20%。由于采樣點(diǎn)每3個(gè)點(diǎn)在同一下風(fēng)向距離，圖5三點(diǎn)線段圖可以看出，同一下風(fēng)向距離，正對(duì)泄漏源的位置濃度最高，兩側(cè)濃度呈對(duì)稱衰減分布。對(duì)稱衰減分布的特征，與上一節(jié)定性分析結(jié)論相一致。

比較實(shí)驗(yàn)測量值、數(shù)值模擬值與擴(kuò)散模型計(jì)算值，實(shí)驗(yàn)測量值與數(shù)值模擬值和擴(kuò)散模型計(jì)算值相近，且普遍偏大。圖5三點(diǎn)線段圖中，實(shí)驗(yàn)測量值存在個(gè)別點(diǎn)(采樣點(diǎn)8)有比較大的偏差，質(zhì)量濃度值比數(shù)值模擬及計(jì)算值小，但與同一下風(fēng)向距離的2個(gè)點(diǎn)數(shù)值相近。

2.2.3原因分析

測量值較大，分析可能有2個(gè)原因：爆炸后迅速蒸發(fā)的推進(jìn)劑氣體殘留導(dǎo)致濃度增加；推進(jìn)劑爆炸燃燒改變了局部溫度，改變了大氣的穩(wěn)定度，也改變了氣體擴(kuò)散能力，使得濃度增加，擴(kuò)散模型計(jì)算中擴(kuò)散參數(shù)應(yīng)該進(jìn)行修正。

圖5三點(diǎn)線段圖中，存在個(gè)別點(diǎn)有比較大的偏差，且在同一下風(fēng)向三點(diǎn)中，峰值減小。分析認(rèn)為，實(shí)際測量過程中風(fēng)向會(huì)有變動(dòng)，局部氣流會(huì)影響到氣體擴(kuò)散濃度，而擴(kuò)散模型計(jì)算和數(shù)值模擬計(jì)算忽略這種情況的理想條件，會(huì)導(dǎo)致存在一定的偏差。另外，考慮到測量的方法比較粗糙，也存在偶然誤差和測量誤差的可能。

總體來看，實(shí)驗(yàn)測量值與數(shù)值模擬值及計(jì)算值的相關(guān)度較大，但實(shí)驗(yàn)值偏大，存在多種因素影響，在爆炸燃燒后泄漏擴(kuò)散情形下，應(yīng)對(duì)計(jì)算公式進(jìn)行參數(shù)修正使其滿足安全和精度要求。

3　結(jié)論

1)氣體擴(kuò)散模型與數(shù)值模擬及實(shí)驗(yàn)結(jié)果基本一致，趨勢相同，數(shù)值整體偏?。徽f明擴(kuò)散模型能較好地反映氣體擴(kuò)散的濃度分布，但存在一定偏差。

2)由于推進(jìn)劑燃燒和氧化反應(yīng)，擴(kuò)散區(qū)域溫度上升，大氣穩(wěn)定度降低，故實(shí)際濃度比理論計(jì)算值大。

[1]叢繼信, 王力, 張光有.液體推進(jìn)劑職業(yè)中毒風(fēng)險(xiǎn)評(píng)價(jià)及防護(hù)對(duì)策研究[J].中國安全生產(chǎn)科學(xué)技術(shù),2012,8(7):40-45.

CONG Jixin,WANG Li,ZHANG Guangyou.Study on risk assessment and protective solutions for occupationalpoisoning of liquid propellant[J].Journal of Safety Science and Technology,2012,8(7):40-45.

[2]任建國, 魯順清. 氣體擴(kuò)散數(shù)學(xué)模型在安全評(píng)價(jià)方面的應(yīng)用[J]. 中國安全科學(xué)學(xué)報(bào),2006,16(3):12-16.

REN Jianguo, LU Shunqing. Application of gas diffusion mathematical model to safety assessment[J].China safety science journal, 2006,16(3):12-16.

[3]陳新華, 武江濤, 佟連捷，等. 液體火箭爆炸后毒氣擴(kuò)散研究[J]. 推進(jìn)技術(shù),1999,20(5):6-10.

CHEN Xinhua, WU Jiangtao, TONG Lianjie，et al.Research on Toxic gas diffusion after liquid propellant rocket explosion[J]. Journal of Propulsion Technology, 1999,20(5):6-10.

[4]蔣維楣, 孫鑒濘, 曹文俊，等. 空氣污染氣象學(xué)教程[M]. 北京:氣象出版社,2004.

[5]陳新華, 莊逢辰. 液體火箭推進(jìn)劑在自然環(huán)境中蒸發(fā)特性[J]. 熱科學(xué)與技術(shù), 2005,4(4):304-308.

CHEN Xinhua, ZHUANG Fengchen. Analysis of liquid rocket propellant evaporation characteristic in natural environment[J], Journal of Thermal Science and Technology, 2005,4(4):304-308.

[6]M.Adon, C.Galy-Lacaux, C.Delon， et al. Dry deposition of nitrogen compounds (NO2, HNO3, NH3) sulfur dioxide and ozone in west and central African ecosystems using the inferential method[J]. Atmos. Chem. Phys., 2009(13):11351-11374.

[7]D.B.Turner. Workbook of atmospheric dispersion estimates: an introduction to dispersion modeling[M]. NYC:CRC Press, 1994.

[8]SAFITRI, X.Gao, M. S. Mannan. Dispersion modeling approach for quantification of methane emission rates from natural gas fugitive leaks detected by infrared imaging technique[J]. Journal of Loss Prevention in the Process Industries, 2011,24(2):113-120.

[9]陳新華, 向四桂. 航天發(fā)射場有毒氣體污染范圍安全性研究[J]. 指揮技術(shù)學(xué)院學(xué)報(bào),2001,12(3):59-62.

CHEN Xinhua, XIANG Sigui. Study of Safety in pollution Range of poisonous gas to space vehicle launching site[J]. Journal of Institute of Command and Technology, 2001,12(3):59-62.

[10]陳新華.液體火箭推進(jìn)劑爆炸毒氣逸散理論與實(shí)驗(yàn)研究[R].北京：裝備指揮技術(shù)學(xué)院,2001.

[11]吳玉劍,潘旭海.障礙物地形條件下重氣泄漏擴(kuò)散實(shí)驗(yàn)的CFD模擬驗(yàn)證[J].中國安全生產(chǎn)科學(xué)技術(shù),2010,6(3):13-17.

WU Yujian，PAN Xuhai. Simulation and verification of CFD on dispersion of heavy gas leakage in obstacle terrain[J].Journal of Safety Science and Technology,2010,6(3):13-17.