煤油貯箱冷氦鼓泡增壓過程數值研究

周芮，程光平，張浩，任楓，王舜浩，張小斌

（1 浙江大學制冷與低溫研究所，浙江杭州310027； 2 上海宇航系統工程研究所，上海201109）

引 言

某新型運載火箭擬采用冷氦增壓技術，在國際上中國是第二個掌握這項技術的國家[1-2]。冷氦增壓技術應用于使用低溫推進劑一類的火箭中，氦氣瓶放置在液氧箱中，低溫氦氣進入煤油箱增壓。利用高壓低溫下，氦氣密度大的特點，一個冷氦氣瓶置于液氧箱中的貯氣量約為同壓力、同容積下常溫氣瓶貯氣量的3倍。因此該技術可以大大減輕貯箱系統質量，提高火箭運載能力。根據新型運載火箭設計方案，擬采用低溫氦氣在液面以下以直接鼓泡方式噴出，然后與煤油換熱升溫后進入氣相空間，提高氣體溫度和增壓能力。但這種煤油貯箱增壓方法存在與冷氦氣流接觸的局部煤油溫度降低結冰的風險，結冰沉底將堵塞發動機濾網，從而引起發動機故障。另外煤油排出過程可能夾雜氦氣泡進入發動機，誘導發動機泵發生氣蝕，使其不能正常工作。而且氣泡在煤油中上升速度慢，可能存在增壓延遲。針對這些可能的問題，本文對冷氦氣在貯箱中心單孔噴入和環向多孔噴入兩種工況（兩種工況冷氦氣流量相同）開展了計算流體力學（CFD）數值分析，分別建立了軸對稱和局部三維幾何模型。為驗證數值模型，計算結果首先與基于原理性實驗貯箱，冷氦氣在常溫水中的鼓泡實驗觀察進行了對比，并開展了煤油的數值計算和分析。

國內研究者對于冷氦增壓系統展開了許多模擬和實驗研究[3-5]。針對系統的研究，曾源華[6]分析了模擬火箭飛行過程中加熱20.4 K低溫氦的實驗結果。邢力超等[7-8]建立了液氫溫區冷氦增壓系統，采用深低溫制冷機配合高壓低溫換熱貯罐獲得20 K、21 MPa 冷氦氣源，同時使用管殼式加熱器模擬火箭上的換熱器性能，實驗系統安全可靠性得到了證實。針對具體的增壓過程，白曉瑞[9]采用AMESim 軟件對液氧增壓系統進行方案設計，并且比較了氦氣、氮氣和空氣的增壓效果，最終結果表明氦氣的增壓效率最高。王承等[10]利用仿真軟件Sinda/Fluint基于集總參數法對冷氦氣瓶的加注過程和系統增壓過程進行建模和計算分析，結果對氣枕壓力和氣瓶熱力狀態之間的關系進行了深入探討。針對增壓系統中的換熱器研究，杜正剛等[11]對各工況下換熱器的換熱性能和壓力損失展開了實驗和模擬研究，結果表明冷氦氣流量增大會導致換熱器出口溫度降低和壓力損失。趙耀中等[12]對比了三種換熱器，結果表明采用低溫氦氣-常溫氮氣的換熱器最好，原理恰當，容易實現且便于操作。

計算流體力學（CFD）方法也被廣泛應用于氦氣增壓低溫推進劑貯箱（液氧和液氫）熱力學過程的研 究[13-19]。NASA Glenn 研 究 中 心 的Kassemi 等[20]基于Fluent 軟件開展了二維軸對稱液氫貯箱自增壓、循環去熱分層過程的CFD 數值研究。結果發現Lee模型的經驗系數對計算結果影響較小。即使Reynolds 數顯示流動是完全湍流范圍，層流模型計算得到的壓力動態變化比基于湍流模型計算結果與實驗結果更一致,原因是湍流有效熱導過大。國外有研究者基于Fluent軟件模擬內壁帶方形和半圓形肋片的圓柱貯箱的自增壓過程，目的是檢驗肋片結構對貯箱壓力動態變化的影響，研究了方形肋片相對間隔以及肋片結構的影響[21]。分析基于2D 軸對稱模型，兩相流基于VOF 模型，湍流采用二方程standard k-ε 模型。結果顯示半圓形肋片對減小壓力增加速率效果比方形肋片更好。國內低溫推進劑國家重點實驗室聯合西安交通大學基于CFD 方法開展了液氧貯箱在加壓加注階段（頂部噴入高溫氧氣）以及發射過程最后級自增壓的研究[22-23]。所采用的兩相數值模型同樣基于2D 軸對稱及VOF 模型，湍流采用二方程standard k-ε 模型。系統研究了貯箱絕熱結構外環境溫度對自增壓的影響，考慮了空間輻射的影響。研究結果顯示，環境溫度變化（253 ～313 K）對液氫增壓幾乎沒有影響。Zhou 等[24]用Fluent 基于VOF 方法，采用層流模型模擬了地面液氫儲罐的自增壓過程，探討了儲罐臨界壓力對蒸發過程的影響。考慮到儲罐實際的結構，模擬過程中理論推導了兩種儲罐表面（包裹泡沫和不包裹泡沫）是否結冰的熱力學過程，并通過UDF 實現并應用于邊界條件。通過與實驗數據的對比，驗證了該模型的準確性，深入探討了液氫的三個蒸發階段并且為低溫儲罐的設計提供了參考。

以上研究內容都是對室溫氦氣進入低溫貯箱氣相空間的增壓過程進行數值研究，低溫氦氣進入液相鼓泡增壓的數值研究報道較少[25]。對比兩種物理過程，前者氣液間滑移速度(slip velocity)可忽略不計，兩相流動計算都是基于體積含量模型（VOF）。該模型為“單流體”模型，即假設氣液兩相流為一種混合流體，同時不考慮相間的相對速度。而后者低溫氦氣高速噴射進液體，氣液間相對速度不能忽略，因此有必要構建Euler-Euler 模型。相比于VOF模型，Euler-Euler 模型對兩相分別列出控制方程，并考慮相間相對速度，但同時也更復雜，對計算資源要求更高。本文基于Euler-Euler 模型對煤油貯箱中冷氦鼓泡增壓進行了數值模擬，模擬結果與水箱的冷氦鼓泡實驗進行了對比。所得結果希望為某新型運載火箭的煤油增壓方案設計提供參考。

1 數值模型

建模的貯箱直徑1.2 m，圓柱段高度1.8 m，兩端為橢球形封頭，直徑比1.6，底部有排出口。圖1 和圖2 分別給出了氦氣單孔鼓泡和環形鼓泡增壓CFD幾何模型及網格，其中單孔鼓泡模型的網格數為49386 個，中心孔半徑為15 mm，在貯箱內部的氦氣輸送管由于管徑相對較小，沒有被建模。環形鼓泡模型的網格數為1428008 個，環向扇形角為30°，每隔14.4°有一個直徑3 mm 的小孔。在兩種鼓泡方式中，氦氣孔到封頭底部之間的距離均為0.575 m。

考慮到低溫氦氣與煤油之間的溫差和速度差，有必要構建更全面復雜的Euler-Euler 兩相模型。同時為減小模型復雜性，本文通過煤油的溫度以及是否有時間結冰兩步來判斷結冰可能性，因此模型中不包括液固相變，相應的控制方程如下。

質量守恒方程

動量守恒方程

式中，k 指氣相（G）或液相（L）。τD,k為氣液相間拖曳力，拖曳力的表達式為

能量守恒方程

圖1 單孔鼓泡增壓方案儲罐網格Fig.1 Mesh scheme of central injection from a core

圖2 環形鼓泡擴散增壓網格Fig.2 Mesh scheme of circumferential injection from three cores

考慮到流場中既有高Re區域（氦氣噴出口），也存在低Re 區域（遠離噴出口的區域），因此選擇SST k-ω 湍流模型。該模型綜合了k-ω 模型在近壁區計算的所有優點，并在湍流黏度定義中考慮了湍流剪切應力輸運過程，適用更廣，可以用于帶逆壓梯度的流動計算和翼型計算等。

計算中所有的物性數據都來自Refprop 9.1。對于氦氣：密度由理想氣體公式計算。壓力和溫度的影響很小，因此為常數cp=5193 J/(kg·K)。熱導率和動力黏度為溫度的函數。本文首先建模了低溫氦氣在水中的鼓泡增壓過程，以與實驗對比。計算中，水的物性都為常數，參考溫度為20 ℃。

煤油物性也都為常數，其中密度為805 kg/m3，比熱容為1967.8 J/(kg·K)，熱導率為0.149 W/(m·K)，動力黏度為0.001 kg/(m·s)。

模擬中，單孔和環形鼓泡中氦氣和煤油總排放量分別為0.04 kg/s 和55.4 kg/s，其中環形鼓泡方案中的每個小孔流量為0.0016 kg/s。整個流體域初始溫度為280 K，絕對壓力為0.32 MPa。

基于Fluent18.1 進行一階隱式非穩態數值計算。迭代方案：壓力-速度耦合方法Coupled，相比于壓力和速度Separated 方案，Coupled 方案同時求解動量方程和質量守恒方程，同步收斂壓力和速度，因此在計算非穩態時能夠更快地收斂。其他方程的離散格式設置如下：壓力求解，PRESTO!；動量方程，二階迎風格式，體積含量，一階迎風格式；能量方程，一階迎風格式。方程收斂標準：質量守恒方程殘差RRR＜10-3，所有速度量殘差RRR＜10-3，能量殘差RRR＜10-6。

2 模型驗證

圖3 冷氦環形鼓泡實驗貯箱內部Fig.3 Interior of ringlike holes scheme tank

圖4 冷氦在水中環形鼓泡實驗快照Fig.4 Experiment moment during cold helium pressurization process in water tank

為驗證方案可行性及數值模型，上海宇航系統工程研究所和上海交通大學黃永華課題組合作開展了90 K 氦氣在水中環形鼓泡實驗，圖3 給出了內部相機從上往下拍攝的實驗裝置內部，包含冷氦進入管、熱電偶支架以及環形排氣管等，進入管和環形管通過三根水平方向的旁管聯通。圖4給出了另一個攝像機從水平方向拍攝的環形排氣管附近冷氦氣剛開始排放時的現象。從圖4 可以看出，環形管上各個小孔流量并不均勻，靠近旁通管的小孔首先噴出氦氣，而后其他小孔噴出氣體。氦氣噴出小孔后，馬上與相鄰小孔噴出的氣體混合，形成大氣團，攜帶液體上升，即存在明顯的氣液混合和攜帶現象。

圖5 給出了以小孔為中心，左右偏轉7.5°三個截面上典型的冷氦氣鼓泡氣液相分布，7.5°的偏轉截面處于兩個小孔的中間。由于三個截面都出現了氦氣，可見冷氦氣噴入水中后擴散并與相鄰小孔的氦氣進行了混合，與實驗現象一致。另外布置的熱電偶溫度傳感器由于離環形管較遠，實驗過程沒有觀察到溫度明顯下降，與數值結果一致。實驗過程所有管路上沒有觀察到結冰現象。

3 計算結果分析

3.1 單孔鼓泡

圖5 環形鼓泡狀態下模擬的不同截面上冷氦氣從水中鼓泡噴出結果（“1”指的是全氣相）Fig.5 Simulation results of cold helium pressurization process in water tank on three different sections (legend denotes volume fraction of gas)

圖6 單孔鼓泡時儲罐中的壓力隨時間的變化曲線Fig.6 History of ullage pressure in one hole scheme tank

圖6給出了中心軸上單孔鼓泡時貯箱氣相空間壓力隨時間的變化。考慮到實驗中以水代替工質，因此分別模擬了液相為水和煤油兩種工況，以定量說明水是否可以代替煤油得到一致的實驗結果。首先，由于排出液體體積流量和氦氣進入流量完全相同，因此兩者的壓力變化趨勢幾乎一致，但數值上水儲罐的壓力始終比煤油儲罐壓力稍小，原因可能是水罐氣相空間中氦氣溫度比煤油儲罐低。圖7給出了t=36 s時刻水儲罐和煤油儲罐氣枕空間的溫度云圖分布，從圖中可以明顯看到水儲罐氣枕空間內，氦氣的平均溫度更低。冷氦氣從噴口出來后，與高溫液體進行強烈混合對流和換熱。由于煤油密度比水小（0.8∶1），氦氣在煤油中的浮力比在水中更小，到達液面的時間更長，計算顯示，當工質為水時，最大速度為5.52 m/s，煤油中最大速度為5.24 m/s，因此煤油和氦氣混合和換熱更充分，導致氦氣到達液面后溫度更高。其他結果，如溫度、相分布等定性一致，因此下面只對氦氣/煤油數值結果進行分析。

圖8 給出了氦氣/煤油單孔鼓泡在t=4,16,32 s時刻的溫度分布云圖。 由圖可見，由于(ρcp)氦氣?(ρcp)煤油，氦氣在與液體換熱后溫度很快就上升到了高溫，自始至終低溫區域被限制于氦氣入口處。氦氣相對于煤油高速噴出，因此是強制對流換熱，同時存在強烈的湍流動量交換，在出口附近被冷卻的煤油被氦氣帶走，同時更遠的熱流體流過來補充，溫度場達到動態平衡。此過程中，由于熱傳導引起的熱擴散相對于對流換熱影響較小，因此，溫度場沒有隨著時間延長擴散。同時，受煤油側的流動和壓力場擾動，以及強烈湍流的影響，溫度場分布并不能保持不變的分布。

圖7 t=36 s時刻水儲罐和煤油儲罐氣枕空間的溫度云圖分布Fig.7 Temperature contours of gas phase inside water/kerosene tank at t=36 s

圖8 氦氣/煤油單孔鼓泡在t=4,16,32 s時刻的溫度分布云圖Fig.8 Temperature contour of helium/kerosene in central injection scheme tank at t=4,16,32 s

圖9 在t=4 s時刻x=0.664 m處的溫度和液相含量分布曲線Fig.9 Temperature and liquid fraction lines at x=0.664 m and t=4 s

圖9給出了在t=4 s時刻在噴口上方0.006 m（即x=0.664 m 處）沿水平方向（如圖8 實線所示）的溫度和相含量分布。要使結冰現象發生需要滿足三個必要條件：①有液體存在；②溫度低于固化溫度（這里煤油固化溫度為T=213.15 K）；③時間足夠長。從圖中可以看出，滿足①和②兩個條件的區域在0.015 m≤y≤0.019 m，寬度約4 mm，該區域對應的液相含量約為0～0.04。因此該區間內的液滴是否結冰決定于第③個條件。當噴口附近的煤油液滴溫度降到固化溫度213.15 K（-60℃），要使其進一步固化，則溫度需要從液滴外表面開始到中心進一步降低（考慮到一定的過冷度）。忽略液滴內部對流換熱，從外表面到球心為一維導熱，則根據傅里葉系數（Fo）

3.2 環形鼓泡

圖10 給出了環形鼓泡時貯箱氣枕壓力隨時間的變化，為對比，圖中同時給出了中心孔鼓泡時的氣枕壓力變化。由圖可見，相比于中心孔進氣，環形進氣時氣枕壓力存在明顯的滯后，且隨著時間增長壓差越來越大。原因是環形進氣時每個小孔流速更小，因此到達氣枕的時間更長，同時存在與液體的混合與旋渦（圖5 和圖11），部分氣體長時間滯留在液相中，所以環形進氣時氣枕壓力上升緩慢。若要讓儲罐達到相同的壓力，在采用環形鼓泡方案時，則需要進一步提高噴氣壓力。

圖11 給出了t=0.3,0.9,2.4,12.5 s 時，中間的排氣孔中心截面的氣相云圖。相應的溫度分布與中心孔鼓泡類似，因此認為不存在結冰的可能性。由圖可以看出排出的氦氣逐漸往上移動，隨時間增加，氦氣在液相中擴散，所占體積越來越大，對液相存在強烈的攪拌作用，這現象與實驗一致。計算結果表明在t=7.5 s 左右氦氣到達氣液界面進入氣相空間。

圖10 環形和中心孔鼓泡氣枕壓力隨時間的變化Fig.10 History of ullage pressure of one hole and ringlike holes schemes

為進一步證實氦氣的擴散過程，圖12 給出了t=0.3,0.6,0.9 s 時環形小孔上方水平平面上氣相體積含量云圖。由圖可清楚地看到氦氣在水平方向上的擴散，大概在0.6 s 時，相鄰的云團開始聚合。一個問題是，為什么中心鼓泡時氦氣沒有擴散，環形鼓泡時出現了嚴重的擴散現象。定性的解釋是中心孔鼓泡時雖然噴出孔面積較大，但氣液之間的接觸面積要小得多，相同體積流量時噴出速度小，對貯箱煤油的流動攪拌局限在小孔上方的圓柱形體積內。但是環形鼓泡時，小孔總面積要小得多，噴出速度相對更快，氣液接觸面積也大得多，氣體更容易被液體攜帶，且相鄰小孔對流場的影響區域存在重疊，加劇了氦氣的擴散。

4 結 論

本文基于Euler-Euler 兩相數值模型和SST k-ω湍流模型，對冷氦氣在煤油貯箱中的兩種鼓泡結構中心孔式和環形小孔式進行了數值模擬。模型與氦氣在水中的環形小孔鼓泡實驗結果進行了定性驗證。對比了煤油和水鼓泡時氣枕壓力的變化。基于數值結果對煤油結冰可能性進行了理論分析。得出主要結論如下。

（1）構建的數值模型定性準確地捕捉到了兩種鼓泡方式的流動和傳熱特性。

（2）數值結果表明低溫氦氣進入液體后在有限行程內馬上升溫到接近煤油溫度，理論分析表明固化可能性微乎其微，實驗觀察證實了這個結論。

（3）中心孔鼓泡中氦氣幾乎沒有擴散，在噴射孔上方至液面下形成半徑比中心孔稍大的氣柱。但是環形小孔鼓泡方案數值和實驗都觀察到了嚴重的氦氣擴散現象，不滿足煤油排出要求。

圖11 t=0.3,0.9,2.4,12.5 s時中間的排氣孔中心截面的氣相云圖Fig.11 Phase contour of middle section of tank at t=0.3,0.9,2.4,12.5 s

圖12 t=0.3,0.6,0.9 s時環形小孔上方水平平面上氣相體積含量云圖Fig.12 Gas fraction contour of horizontal plane above ringlike holes at t=0.3,0.6,0.9 s

（4）由于氦氣和煤油換熱更充分，相同條件下煤油儲罐氣枕壓力比水儲罐氣枕壓力上升更快，但差別不大，可以用水代替煤油進行冷氦氣鼓泡實驗。