基于VOF模型的偏二甲肼泄漏液池擴展過程的數值模擬

詹 翔，崔村燕，周宵燈，趙繼廣，辛騰達

(1.航天工程大學研究生管理大隊，北京 101416；2.航天工程大學宇航科學與技術系，

北京 101416；3.航天工程大學電子與光學工程系，北京 101416)

偏二甲肼(Unsymmetrical Dimethylhydrazine,UDMH)易燃、易爆、劇毒、腐蝕性強，在航天發射場使用的安全問題歷來受到人們的高度關注[1]。肼類燃料對非金屬的溶脹作用很強烈，因此法蘭連接處的密封、閥門和泵填函部位的密封長時間使用后容易失效[2]。實際上在偏二甲肼儲存和加注期間，容器閥門法蘭腐蝕、設備失效、焊接缺陷以及人員的誤操作等原因[3]都會造成推進劑不同程度的泄漏。偏二甲肼泄漏后在地面上形成液池，遇到火星易著火爆炸，也容易蒸發形成有毒氣體接觸到人員和設備,可能對其造成損害[4]，嚴重威脅著周圍設備和工作人員的安全。而研究偏二甲肼泄漏液池的擴展規律是進一步研究推進劑蒸氣擴散規律、開展推進劑安全評價工作、制定推進劑泄漏事故應急預案等的重要前提。在以往對偏二甲肼泄漏液池擴展過程的研究中[5]，對偏二甲肼泄漏液池做了極簡處理，即簡單地以圓形液池代替，這樣并不能真實地反映偏二甲肼泄漏液池的分布規律。為此，本文采用VOF模型對偏二甲肼泄漏形成液池的擴展過程進行了數值模擬，研究了在受限空間內外泄偏二甲肼液池的擴展規律。

1 數值模擬方法

1. 1 VOF模型

數值模擬軟件FLUENT中的多相流模型可以模擬不同流體間的相互流動，其中VOF模型通過自由液面的重構來獲得自由液面的位置和形狀[6]，是處理液氣兩相交界、互相推動力瞬變問題的重要方法[7]。VOF模型采用k-ε紊流模型，其基本的控制方程包括連續方程、動量方程、湍流方程和湍動能耗散率方程[8]。利用VOF模型追蹤液體的自由表面是數值模擬方法中較為有效的方法，目前已有大量的研究，如陳群等[9]通過對階梯溢流壩建模，對溢流壩面流場的紊流進行了數值模擬，得到的數值模擬結果與實測值相比分布規律近似，數據相差較小；關大瑋等[10]采用FLOW3D軟件，將RNGk-ε模型和VOF模型結合，對溢流壩三維水流進行了數值模擬，模擬結果對照縮比物理模型試驗結果大體上吻合；史曉蒙等[11]通過建立庫區地面輸油管道泄漏擴散過程的三維CFD模型，采用FLUENT軟件中的VOF模型對油品泄漏擴散過程進行了數值模擬，分析了油品泄漏后的擴散趨勢，該研究結果有助于分析泄漏事故后果影響范圍以及合理制定防護及救援措施。

1. 2 自由液面的重構

利用數值模擬方法研究偏二甲肼泄漏液池的擴展規律，重點在于處理偏二甲肼和空氣的相交界面。FLUENT軟件中的控制體積公式需要計算穿過控制體積面的對流和擴散通量，并且與控制體積內部的源項平衡[12]。劃分好網格后FLUENT軟件即可按給定的初始運動界面開始計算，但是能捕捉到的新運動界面是一個模糊的輪廓。為了提高界面的分辨率和精度，FLUENT軟件必須在每個時間步長計算完成后對界面進行重構。對于VOF模型，FLUENT軟件中有4種方案計算流體界面的通量：幾何重構方案、施主受主方案、歐拉顯式方案和隱式插補方案[13]。其中，幾何重構方案求解簡單，本文主要采用幾何重構方案。

在幾何重建方案中，若單元被某一項填滿，那么FLUEN軟件將以標準插值方案計算流體界面的通量，如果單元包含兩相之間的交界面，則使用幾何重建方案。幾何重建方案利用分段線性的方法描繪流體之間的界面，見圖1。

圖1 真實流體界面與幾何重構流體界面Fig.1 Real fluid interface & geometric reconstruction

由圖1可以明顯地看出，真實流體的界面是圓滑的弧形界面[見圖1(a)]；而幾何重構方案則將原來的圓滑弧形界面在單元格內由線段代替[見圖1(b)]。

假定兩種流體之間的界面在每個單元內有線性的斜面，并使用這個線性界面形狀來計算穿過單元面之間的對流，在計算流體界面通量時，首先根據體積分數和單元的相關信息來計算出線性界面的位置，再通過界面上的切向速度分布信息來計算每個流體界面的通量。

當實際交界面的尺度相對于網格尺度很小時，幾何重構出的線性界面將無法反映真實的流體交界面，見圖2。

圖2 交界面與網格尺寸的對比圖Fig.2 Comparison of interface and grid size

由圖2可見，氣泡的大小比單個計算單元格還小，若按照幾何重構方案，在一個單元格內用一條線段代替交界面，則無法正常刻畫氣泡的輪廓。Rabha等[14]認為網格數與氣泡直徑的比值至少為16時才能較好地刻畫出氣泡的輪廓。本文為了較好地刻畫偏二甲肼泄漏過程中出現的液滴，在局部區域做了網絡加密處理。

2 幾何模型及邊界條件設置

以某偏二甲肼容器參數為依據，對容器實際外形條件適當簡化后建立幾何模型。某偏二甲肼容器高35 m、半徑為1.675 m，假定加泄連接閥處法蘭腐蝕發生連續泄漏，并通過容器外壁面上2 m處的連接口向外泄漏。由于泄漏口以上的部分在實際計算中對偏二甲肼的泄漏幾乎沒有影響,因此在不影響數值模擬結果的前提下對幾何空間進行了合理簡化，只留下偏二甲肼容器的部分弧形壁面和泄漏口正前端的部分空間，見圖3。

圖3 總體結構與計算區域Fig.3 Structure and calculation area

根據實際要求，受限空間環境參數需控制在溫度為15℃，相對濕度為50%，風速為0.5 m/s。以偏二甲肼容器加泄連接閥處法蘭接觸面發生泄漏為例，偏二甲肼泄漏量的計算公式為[15]

(1)

式中:Q為偏二甲肼泄漏量(cm3/s)；h為間隙的大小(cm)；Δp為設備內外的壓力差(MPa)；η為偏二甲肼的黏度(Pa·s)；r1為法蘭內徑(cm)；r2為法蘭外徑(cm)。

偏二甲肼容器內部為正常大氣壓，求得其泄漏速率為45.8 cm3/s。假定忽略泄漏過程中容器內壓力的變化，則泄漏速率保持不變，進而可根據加泄連接口面積大小獲得出口速度為0.03 m/s。本次模擬假設偏二甲肼容器壁面為絕熱、無滑移，且流體區域不可壓縮，呈湍流狀態，湍流分布符合伯努利假設；由于空間有限，假設空間內空氣分布均勻，不考慮空氣的分層和穩定度。

3 偏二甲肼泄漏液池擴展過程的數值模擬與分析

3. 1 偏二甲肼泄漏液池的分布

不同時刻偏二甲肼泄漏液池分布的數值模擬結果見圖4。

圖4 不同時刻偏二甲肼泄漏液池的分布Fig.4 UDMH liquid pool distribution at different times

由圖4可見，當t=0.15 s時，偏二甲肼在泄漏口附近聚集,在重力作用下，一部分偏二甲肼從泄漏液主體分離，單獨沿壁面下落;當t=1.50 s時，下落的偏二甲肼泄漏液分裂成面積約為3 cm2的液滴，散布在壁面上；當t=3.00 s時，偏二甲肼泄漏液的一部分已經到達地面；當t=4.05 s時，偏二甲肼在地面上聚集形成面積大小約為224.35 cm2的液池。

3. 2 垂直壁面上偏二甲肼泄漏的流態分析

圖5給出了不同時刻垂直壁面上偏二甲肼泄漏液池面積的增長曲線圖。

圖5 不同時刻垂直壁面上偏二甲肼泄漏液池面積的增長曲線圖Fig.5 Growth curve of the UDMH liquid pool area on the vertical wall at different times

由圖5可見，在垂直壁面上偏二甲肼泄漏液池的擴展過程大致可分為A、B、C三個過程：A過程中，偏二甲肼泄漏的時間較短，泄漏量較少，泄漏液體積小，液體黏性力和表面張力對液池擴展起到了一定的阻礙作用[16]，此時泄漏液受重力的加速作用還不明顯，在泄漏口附近聚集并緩慢下落，液池面積增長速度比較慢；B過程中，一部分液體在重力的作用下獲得了相對高的下落速度，在壁面摩擦力的作用下，分散成多個小液滴，小液滴在壁面上橫向分散，因此該過程中液池擴展速度迅速增加，且在B過程末期，偏二甲肼泄漏液不再形成新的液滴，壁面上偏二甲肼液池呈現下寬上窄的分布規律，此時液池的面積增長速率也逐漸放緩；偏二甲肼泄漏液到達地面時，開始進入C過程，該過程中由于偏二甲肼液滴到達地面且無新的液滴產生，垂直壁面上偏二甲肼泄漏液池的總面積不斷減少，且C過程結束后所有分散的液滴達到地面，當垂直壁面上只剩下垂直向下且細長的偏二甲肼泄漏液池時，偏二甲肼泄漏液池的面積將基本不再變化。

圖6給出了偏二甲肼泄漏液池具有A、B、C三個過程典型特征的流態圖。

圖6 不同時刻垂直壁面上偏二甲肼泄漏的流態圖Fig.6 Flow pattern of the UDMH leakage on the vertical wall at different times

當垂直壁面上偏二甲肼泄漏液池的面積不再變化后，由數值模擬結果可得到偏二甲肼泄漏液池的面積約為301 cm2。偏二甲肼蒸發速率的計算公式為[17]：

vu=2.603 6Aγ0.553×(0.001 935T2+0.018 988T+0.511 81)×(-2H2+4.46H+0.026 5)

(2)

式中：vu為偏二甲肼的蒸發速率(mg/min)；A為偏二甲肼泄漏液池的面積(m2)；γ為風速(m/s)；T為溫度(℃)；H為相對濕度。

由以上環境條件，可計算出在垂直壁面上偏二甲肼的蒸發速率約為0.119 mg/min。

3.3 偏二甲肼泄漏在地面上形成的液池的擴展過程分析

偏二甲肼泄漏到達地面后緊靠容器壁面形成半圓形的液池，圖7給出了不同時刻地面上偏二甲肼泄漏液池面積的增長曲線。

圖7 不同時刻地面上偏二甲肼泄漏液池面積的增長曲線圖Fig.7 Growth curve of the liquid pool area of UDMH leakage on the ground

由圖7可見，最初少量接觸到地面的偏二甲肼泄漏液匯聚在一起，由于液體自身黏性力和地面摩擦力的阻礙作用,此時偏二甲肼泄漏液池面積的增長速度較慢；隨著偏二甲肼在地面上聚集量的增加，偏二甲肼泄漏液池在自身重力和下落滴液沖擊力的作用下克服液體自身黏性力和地面摩擦力，形成面積更小的圓形液池向四周擴散，此時偏二甲肼泄漏液池面積的增長速度也有所增加。圖8給出了不同時刻地面上偏二甲肼泄漏液池的擴展過程。

圖8 不同時刻地面上偏二甲肼泄漏液池的擴展過程Fig.8 Expansion process of the liquid pool of UDMH leakage on the ground at different times

3. 4 泄漏速率對偏二甲肼漏泄漏液池擴展的影響

為了研究不同泄漏速率對偏二甲肼泄漏液池擴展規律的影響，將泄漏速度加倍即出口速度改為0.06 m/s，其他條件不變進行了數值模擬計算，得到同時刻(3.00 s)不同泄漏速度條件下垂直壁面上偏二甲肼泄漏液池的擴展情況對比，見圖9。

圖9 同時刻不同泄漏速度下垂直壁面上偏二甲肼 泄漏液池的擴散情況對比Fig.9 Comparison of the expansion of UDMH liquid pool on container wall at the same time with different leakage velocity

由圖9可見，不同泄漏速度對垂直壁面上偏二甲肼泄漏液池的形狀影響較小，在泄漏速度較大的情況下，垂直壁面上分散的小液滴數量明顯減少。由此可見，不同的泄漏速度條件下，垂直壁面上仍然是形成垂直向下且細長的偏二甲肼泄漏液池。

圖10為不同泄漏速度條件下地面上偏二甲肼泄漏液池的擴展面積對比。

由圖10可見，泄漏速度對地面上偏二甲肼泄漏液池擴展面積的影響十分明顯： 當泄漏速度為原來的兩倍時，偏二甲肼泄漏液池的擴展面積是原速度情況下的4倍左右；當泄漏量較大時，偏二甲肼泄漏液池更加接近半圓形，并且液池向四周擴散時分裂形成小面積液池的過程進行得更快。這主要是因為泄漏量的增加，導致泄漏液能獲得更大的動能，以克服地面摩擦力和液體自身黏性力向四周擴散。

圖10 不同泄漏速度下地面上偏二甲肼泄漏液池擴展面積的對比Fig.10 Comparison of the expansion area of UDMH liquid pool on the ground with different leakage velocity

4 結 論

本文通過對偏二甲肼泄漏液池的擴展規律進行數值模擬與分析，得到如下結論：

(1) 偏二甲肼泄漏時，在垂直壁面上主要形成垂直向下且細長的液池；液池前端部分在壁面摩擦力的作用下分散成小的液滴，由于液滴向左右擴散，垂直壁面上偏二甲肼泄漏液池呈現下寬上窄的分布規律；當垂直壁面上所有液滴均到達地面后，壁面上的液池垂直向下流動且寬度基本一致，其面積約為301 cm2。根據環境條件，可求出該液池的偏二甲肼蒸發速率為0.119 mg/min。由此可知，在對偏二甲肼泄漏進行危害性分析時，不能忽視垂直壁面上液池的蒸發作用。

(2) 偏二甲肼泄漏時，地面上的泄漏液在偏二甲肼容器壁面周圍聚集，形成近似半圓形的液池；隨著泄漏液的不斷積累，半圓形液池向周圍擴展的過程中會分裂出面積較小的圓形液池。

(3) 偏二甲肼泄漏時，當泄漏速度加大后，垂直壁面上零散的液滴明顯減少，地面上形成的液池更加接近規則的半圓形;當泄漏速度為原來的兩倍時，地面上偏二甲肼泄漏液池的擴展面積是原速度情況下的4倍左右。