微重力下低溫貯箱內(nèi)推進(jìn)劑相變仿真模型研究

王 夕，王 玨，容 易，黃 輝

?

微重力下低溫貯箱內(nèi)推進(jìn)劑相變仿真模型研究

王 夕1，王 玨2，容 易1，黃 輝1

（1. 北京宇航系統(tǒng)工程研究所，北京，100076；2.中國(guó)運(yùn)載火箭技術(shù)研究院，北京，100076）

低溫推進(jìn)劑具有沸點(diǎn)低、易汽化的特點(diǎn)，相變是低溫推進(jìn)劑長(zhǎng)時(shí)間在軌蒸發(fā)量控制問(wèn)題中需要考慮的首要影響因素。相變模型對(duì)低溫推進(jìn)劑蒸發(fā)仿真起到重要作用，構(gòu)建合理的相變模型成為低溫流體蒸發(fā)量仿真重要的研究方向?；?種相變理論，采用FLUENT軟件二次開(kāi)發(fā)的方法，建立基于相平衡和非平衡理論的4種相變模型，開(kāi)展微重力下液氫推進(jìn)劑蒸發(fā)的數(shù)值模擬，并與國(guó)外探空火箭試驗(yàn)進(jìn)行比較和驗(yàn)證。研究結(jié)果表明：比較4種相變模型對(duì)貯箱內(nèi)壓力升高速率預(yù)示的準(zhǔn)確性，得出了適用于微重力下低溫推進(jìn)劑仿真的相變模型。

低溫推進(jìn)劑；相變；微重力

0 引 言

在液體運(yùn)載火箭的推進(jìn)劑中，低溫推進(jìn)劑具有比沖高、沸點(diǎn)低、難于貯存的特點(diǎn)，因此限制了其長(zhǎng)時(shí)間的在軌使用。對(duì)低溫推進(jìn)劑長(zhǎng)時(shí)間在軌蒸發(fā)量控制問(wèn)題的研究，可以采用試驗(yàn)和仿真的方法。鑒于低溫和微重力試驗(yàn)的難度，基于已有商業(yè)軟件的數(shù)值仿真是對(duì)貯箱內(nèi)低溫推進(jìn)劑開(kāi)展研究的有效途徑。

目前，低溫推進(jìn)劑蒸發(fā)的數(shù)值研究受到廣泛的關(guān)注。在低溫推進(jìn)劑蒸發(fā)量的所有影響因素中，相變是重要的考慮因素，相變過(guò)程的建模是蒸發(fā)量控制研究的焦點(diǎn)問(wèn)題，而發(fā)展相變數(shù)值模型的基礎(chǔ)是描述相變物理過(guò)程的理論模型。相變的理論模型包括相平衡模型和非平衡模型。相平衡模型中，假設(shè)液相和汽相之間達(dá)到平衡時(shí)，界面各相之間的溫度相等，處于平衡狀態(tài)[1]；而在非平衡模型中，Schrage[2]使用基于Maxwell速度分布來(lái)分析液/汽界面的傳質(zhì)過(guò)程，將界面的相變考慮為兩相界面對(duì)液/汽分子的捕獲和逃逸的非平衡過(guò)程；李震東[3]、Tanasawa[4]等總結(jié)了發(fā)生在液/汽相界面處的相變傳質(zhì)理論，歸納為相平衡模型和非平衡模型；劉秋生[5]等使用簡(jiǎn)化的Hertz-Knudsen方程對(duì)熱毛細(xì)作用主導(dǎo)的流動(dòng)進(jìn)行了研究，對(duì)微重力下的相變問(wèn)題具備參考價(jià)值。

基于以上兩類相變理論，建立了相應(yīng)的仿真模型，并用于低溫貯箱推進(jìn)劑蒸發(fā)的仿真中。其中，相平衡模型在低溫貯箱CFD仿真中應(yīng)用廣泛。程向華[6]等根據(jù)相平衡建立了液氧貯箱的二維模型，分析了液氧熱分層的形成過(guò)程及原因；Zilliac[7]等對(duì)貯箱進(jìn)行建模，根據(jù)熱力學(xué)平衡原理建立相變模型。相比于相平衡模型在低溫推進(jìn)劑相變仿真的廣泛使用，非平衡模型的使用很少。Olga[8]等應(yīng)用非平衡模型開(kāi)展了液氫貯箱內(nèi)徑向噴霧的仿真，但并沒(méi)有將非平衡模型與相平衡模型進(jìn)行對(duì)比研究。

實(shí)際的相變物理過(guò)程在界面處存在溫度跳變和非平衡的瞬態(tài)作用，與相平衡模型存在一定的差異，這些差異對(duì)低溫推進(jìn)劑蒸發(fā)數(shù)值模擬仿真準(zhǔn)確性的影響仍未被討論，比較不同相變仿真模型的準(zhǔn)確性和適用性對(duì)低溫流體仿真具有重要作用。

本文基于已有的相變理論展開(kāi)了CFD建模工作，通過(guò)FLUENT用戶自定義函數(shù)，建立基于不同相變理論的蒸發(fā)模型，并通過(guò)適當(dāng)?shù)募僭O(shè)修正蒸發(fā)模型，開(kāi)展了貯箱內(nèi)液氫蒸發(fā)仿真，與國(guó)外Areobee探空火箭試驗(yàn)[9]進(jìn)行了對(duì)比，并分析了不同的相變模型在低溫流體仿真中的應(yīng)用。

1 數(shù)值模型

本文參照Aerobee探空火箭液氫蒸發(fā)試驗(yàn)[9]開(kāi)展了數(shù)值模擬的對(duì)比研究，基于相平衡理論和3種非平衡理論，通過(guò)界面假設(shè)和蒸發(fā)假設(shè)，建立了仿真軟件FLUENT下的相變模型，模擬試驗(yàn)貯箱內(nèi)壓力上升，驗(yàn)證并比較不同模型的結(jié)果，得到適用于微重力下低溫推進(jìn)劑蒸發(fā)量預(yù)示的數(shù)值仿真模型。

1.1 控制方程

采用二維軸對(duì)稱模型，對(duì)液氫貯箱推進(jìn)劑受熱蒸發(fā)的過(guò)程進(jìn)行分析。控制方程為連續(xù)方程、N-S方程和能量方程，對(duì)于每一相，有：

式中為密度；為速度；為動(dòng)力粘度；為時(shí)間；為壓力；為能量；為系統(tǒng)加速度；為溫度；為體積力；為能量源項(xiàng)。

液相密度采用Boussinesq近似，其余物性作為溫度的函數(shù)分段插值，氣相則采用理想氣體模型。

（4）

考慮到蒸發(fā)過(guò)程流速小，因此流動(dòng)采用層流模型。

微重力下需要考慮表面張力作用，采用連續(xù)表面張力模型如下：

1.2 相變模型

相變模型是本文主要研究對(duì)象，基于兩種假設(shè)，開(kāi)展對(duì)4種相變模型的研究。

a）界面假設(shè)：相變僅發(fā)生在相界面處；

b）蒸發(fā)假設(shè)：只有在液相溫度高于飽和溫度時(shí)，才發(fā)生相變，不考慮冷凝的發(fā)生。

1.2.1 模型1

考慮使用相平衡模型[6]，液相的超過(guò)飽和溫度的能量將全部轉(zhuǎn)化為相變，模型1方程為

1.2.2 模型2

Schrage[2]推薦的非平衡模型，即Hertz-Knudsen方程為

凈質(zhì)量流率為正值，表明氣液界面上有凈的質(zhì)量從氣相變?yōu)橐合?，即相界面上發(fā)生凝結(jié)；反之，如果凈質(zhì)量流率為負(fù)值，則相界面上將發(fā)生蒸發(fā)。

1.2.3 模型3

劉秋生[5]等提出采用非平衡熱力學(xué)方法計(jì)算表面張力主導(dǎo)流動(dòng)下的界面蒸發(fā)流量。蒸發(fā)界面假設(shè)為不變形，采用線性近似后的Hertz-Knudsen方程來(lái)描述蒸發(fā)界面的蒸發(fā)流量，即模型3：

1.2.4 模型4

采用相變模型[4,10]，即模型4為

根據(jù)Tanasawa[4]的推薦值及實(shí)際計(jì)算的結(jié)果比較，取系數(shù)為0.2。

1.3 網(wǎng)格無(wú)關(guān)性驗(yàn)證

對(duì)采用的二維軸對(duì)稱網(wǎng)格（見(jiàn)圖1）開(kāi)展網(wǎng)格無(wú)關(guān)性驗(yàn)證（見(jiàn)圖2），使用相變模型1，比較3種不同網(wǎng)格數(shù)量，網(wǎng)格數(shù)分別為5000，10 000和20 000。

圖1 采用的網(wǎng)格

圖2 網(wǎng)格無(wú)關(guān)性驗(yàn)證曲線

3種網(wǎng)格壓力上升的偏差較小，其中與較高網(wǎng)格數(shù)量相比，中等網(wǎng)格數(shù)量的最終壓力差值占總壓力升高的4.17%，該偏差可以被忽略，結(jié)果較為一致?？紤]到網(wǎng)格數(shù)的增加對(duì)計(jì)算效率的影響，采用中等網(wǎng)格數(shù)量開(kāi)展數(shù)值仿真。

2 結(jié)果與分析

對(duì)4種相變模型開(kāi)展研究，對(duì)比不同模型下貯箱內(nèi)壓力上升速率與試驗(yàn)實(shí)測(cè)值的偏差，同時(shí)對(duì)貯箱內(nèi)的流體行為開(kāi)展研究。仿真條件如表1所示，選取與試驗(yàn)相同的初始?jí)毫?，并將初始溫度設(shè)置為該壓力下飽和溫度。假設(shè)貯箱為均勻受熱，選取試驗(yàn)中平均熱流密度作為壁面熱流邊界條件。

表1 數(shù)值模擬條件

Tab.1 Conditions of Simulation

性能參數(shù) 球形貯罐直徑/cm22.86 體積充填率34.3% 平均熱流密度/（W·m-2）473.19 重力/（m·s-2）0.01 初始?jí)毫?MPa0.1241 總受熱時(shí)間/s237

由于采用了兩種假設(shè)，本文將對(duì)界面假設(shè)和蒸發(fā)假設(shè)的合理性展開(kāi)討論，并基于兩種假設(shè)，對(duì)采用的4種相變模型的對(duì)比分析。

2.1 界面假設(shè)

界面假設(shè)將相變界定在界面處，是否采用界面假設(shè)將對(duì)結(jié)果產(chǎn)生影響。采用/不采用界面假設(shè)代表兩種不同的相變?cè)怼?/p>

汽化過(guò)程可分為界面蒸發(fā)和沸騰蒸發(fā)兩種方式。當(dāng)模型不采用界面假設(shè)，液相溫度只要溫度超過(guò)飽和溫度，在任何位置均可發(fā)生汽化，這樣的模型更接近于沸騰。然而，在實(shí)際過(guò)程中，液體內(nèi)部發(fā)生沸騰并產(chǎn)生氣泡需要達(dá)到一定的活化能，并非達(dá)到飽和溫度即開(kāi)始汽化，在較低熱流密度的條件下相變往往以界面蒸發(fā)為主。熱量傳遞到界面需要一定的時(shí)間，當(dāng)界面蒸發(fā)主導(dǎo)相變時(shí)，不采用界面假設(shè)將高估汽化傳質(zhì)速率，即使與實(shí)際的沸騰過(guò)程相比，由于未考慮氣泡生成的過(guò)熱度，不采用加密假設(shè)的模型同樣會(huì)高估汽化傳質(zhì)速率。

采用模型1和模型4來(lái)比較有/無(wú)界面假設(shè)的貯箱壓力上升結(jié)果如表2所示。在237 s內(nèi)，無(wú)界面假設(shè)的模型得到的貯箱壓力和平均壓力升高速率均高于試驗(yàn)值，而采用界面假設(shè)的模型得到的結(jié)果更接近于試驗(yàn)值。結(jié)果表明，不采用界面假設(shè)的模型將高估汽化傳質(zhì)速率，使壓力快速升高；采用界面假設(shè)的模型更為合理，進(jìn)一步支撐了假設(shè)的合理性。對(duì)照試驗(yàn)以界面蒸發(fā)過(guò)程為主，與沸騰過(guò)程有本質(zhì)的區(qū)別，汽化過(guò)程僅發(fā)生在界面處，采用界面假設(shè)的模型也更符合實(shí)際。

表2 界面假設(shè)壓力上升結(jié)果對(duì)比

Tab.2 Pressure Rise Result of Interface Assumption

模型終壓/MPa壓力升高速率/(Pa·s-1) 試驗(yàn)0.7302752557.7 模型1（無(wú)界面假設(shè)）0.8424493031.0 模型4（無(wú)界面假設(shè)）0.8882393224.2 模型1（界面假設(shè)）0.7569942670.4 模型4（界面假設(shè)）0.7682602718.0

在低溫推進(jìn)劑貯箱的相變仿真中，是否采用界面假設(shè)，取決于流體的形態(tài)。若貯箱受熱較小，主要形態(tài)為界面蒸發(fā)，則使用界面假設(shè)；若貯箱受熱較大，主要形態(tài)為沸騰過(guò)程，則應(yīng)當(dāng)考慮不采用界面假設(shè)。

2.2 蒸發(fā)假設(shè)

蒸發(fā)假設(shè)的實(shí)質(zhì)是忽略冷凝的作用。在實(shí)際過(guò)程中，發(fā)生在氣相內(nèi)部（非界面處）的均勻冷凝一般很難發(fā)生，而非均勻冷凝則發(fā)生在有過(guò)冷界面存在的情形下[2]。而試驗(yàn)中不存在大過(guò)冷度的過(guò)冷界面，若認(rèn)為只要低于飽和溫度即發(fā)生冷凝，則可能高估冷凝的作用，造成計(jì)算偏差。因此，不采用純蒸發(fā)假設(shè)可能高估冷凝的作用。

采用模型1和模型3來(lái)比較有/無(wú)蒸發(fā)假設(shè)的貯箱壓力上升結(jié)果如圖3所示。

圖3 蒸發(fā)假設(shè)壓力上升結(jié)果對(duì)比曲線

圖3中，最終的壓力值和平均壓力上升速率可用于比較分析。由于試驗(yàn)中加熱器啟動(dòng)需要一定的時(shí)間，溫度由開(kāi)始加熱到逐漸穩(wěn)定，因此在試驗(yàn)初始階段，加熱的熱流更小。而仿真中，熱流邊界條件設(shè)置為平均熱流密度，并保持恒定值，因而圖中曲線在100 s處壓力比數(shù)值模擬結(jié)果低是合理的，最終壓力是較為合理的比較參數(shù)。

由圖3可知，不采用蒸發(fā)假設(shè)時(shí)，模型3的最終壓力結(jié)果遠(yuǎn)低于試驗(yàn)值，而模型1的最終壓力同樣低于試驗(yàn)值。不采用純蒸發(fā)假設(shè)的模型，高估了冷凝的作用，造成凈汽化速率的降低，壓力上升速率減緩，導(dǎo)致了與試驗(yàn)的偏離。

采用蒸發(fā)假設(shè)的模型1和模型3的最終壓力值與試驗(yàn)值更接近，結(jié)果更準(zhǔn)確。通過(guò)有/無(wú)蒸發(fā)假設(shè)的對(duì)比分析，在沒(méi)有過(guò)冷壁面存在的條件下，采用蒸發(fā)假設(shè)的相變傳質(zhì)模型更合理。

2.3 4種相變模型的結(jié)果分析

在界面假設(shè)和蒸發(fā)假設(shè)下，對(duì)4種蒸發(fā)模型開(kāi)展了對(duì)比研究，以分析模型的準(zhǔn)確性以及相平衡模型和非平衡模型的差別。

4種模型結(jié)果對(duì)比如表3所示，4種模型的壓力結(jié)果如圖4所示。

表3 4種模型結(jié)果對(duì)比

Tab.3 Results of Four Models

模型最終壓力/MPa平均壓力升高速率/（Pa·s-1）偏差 試驗(yàn)0.730 2752557.7— 模型10.756 9942670.44.4% 模型20.657 5342250.8-12% 模型30.734 5462575.70.7% 模型40.768 2602718.06.3%

圖4 4種模型壓力結(jié)果

由表3可知，模型2的平均壓力升高速率偏差最大，為-12.0%。由于在相變模型中考慮了壓力，而壓力是較為敏感的參數(shù)，同時(shí)考慮壓力和溫度也增加了模型的復(fù)雜性，其計(jì)算結(jié)果可能出現(xiàn)較大的偏差。

模型4的壓力升高速率結(jié)果偏差為6.3%，比模型2的偏差更小，比其他兩種模型略為偏大，略高于試驗(yàn)貯箱壓力結(jié)果。模型1是相平衡模型，而其他3種模型均為非平衡模型。模型1的偏差為4.4%，結(jié)果優(yōu)于模型2、模型4，略高于試驗(yàn)貯箱壓力結(jié)果。

模型3的壓力升高速率與試驗(yàn)的偏差最小，偏差為0.7%，最適合用于微重力下液氫貯箱蒸發(fā)量的預(yù)示。

通過(guò)4種模型的對(duì)比和分析發(fā)現(xiàn)，模型1、模型3及模型4的偏差均小于7%，對(duì)微重力下低溫推進(jìn)劑貯箱蒸發(fā)量預(yù)示的偏差均可被接受。其中，模型3考慮了在熱毛細(xì)力主導(dǎo)下的界面蒸發(fā)，對(duì)微重力下的情形更為適用，其結(jié)果相比于其他3種模型具備明顯的優(yōu)勢(shì)，偏差最小，適用于微重力下低溫推進(jìn)劑貯箱蒸發(fā)量預(yù)示和對(duì)貯箱壓力升高的數(shù)值模擬，具備良好的準(zhǔn)確性。

在國(guó)內(nèi)外對(duì)低溫推進(jìn)劑蒸發(fā)的研究中，使用非平衡相變模型的較少，一般使用相平衡模型對(duì)蒸發(fā)過(guò)程進(jìn)行數(shù)值模擬。根據(jù)以上的分析，相平衡模型具備在趨勢(shì)上預(yù)測(cè)貯箱壓力的能力，其仿真結(jié)果偏差在可接受的范圍內(nèi)，這些研究的結(jié)果并不會(huì)在趨勢(shì)上受蒸發(fā)模型的影響。

在微重力這樣的特殊環(huán)境中，模型3相比于相平衡模型和其他非平衡模型具備準(zhǔn)確性優(yōu)勢(shì)。因此，在特定環(huán)境下，采用與環(huán)境相適應(yīng)的非平衡模型將進(jìn)一步提高仿真的準(zhǔn)確性。

3 結(jié) 論

本文利用國(guó)外Aerobee探空火箭試驗(yàn)數(shù)據(jù)對(duì)微重力下液氫貯箱內(nèi)相變仿真開(kāi)展CFD建模工作，根據(jù)不同的相變理論，通過(guò)兩種假設(shè)，對(duì)不同相變模型開(kāi)展比較分析，得到如下結(jié)論：

a）在較低熱流、貯箱內(nèi)以界面蒸發(fā)為主導(dǎo)的情況下，宜采用界面假設(shè)；在沒(méi)有過(guò)冷壁面存在的條件下，對(duì)低溫推進(jìn)劑受熱蒸發(fā)相變的數(shù)值模擬宜采用蒸發(fā)假設(shè)。

b）本文比較了4種相變模型對(duì)封閉液氫貯箱內(nèi)壓力升高的預(yù)測(cè)，其中包括1種相平衡模型和3種非平衡模型，結(jié)果表明，模型1、3、4均能在趨勢(shì)上對(duì)微重力下液氫封閉貯箱壓力升高趨勢(shì)進(jìn)行預(yù)測(cè)，可有效模擬貯箱中的流體形態(tài)和溫度分布。相平衡模型1具備對(duì)貯箱內(nèi)壓力升高趨勢(shì)的預(yù)測(cè)能力。

c）非平衡模型3仿真結(jié)果更接近試驗(yàn)，對(duì)微重力下液氫蒸發(fā)的數(shù)值模擬具備良好的準(zhǔn)確性，在微重力環(huán)境下的仿真中具備優(yōu)勢(shì)。在特定環(huán)境下，采用與環(huán)境相適應(yīng)的非平衡模型將進(jìn)一步提高仿真的準(zhǔn)確性。

[1] 吳雙應(yīng), 曾丹苓, 黃雙, 等. 微層蒸發(fā)模型的非平衡熱力學(xué)分析[J]. 熱科學(xué)與技術(shù), 2005, 4(4): 283-288.

Wu Shuangying, Zeng Danling, Huang Shuang, et al. Non-equilibrium thermodynamic analysis of micro-layer evaporation model[J]. Journal of Thermal Science and Technology, 2005, 4(4): 283-288.

[2] Schrage R W. A theoretical study of interphase mass transfer[M]. NewYork: Columbia University Press, 1953.

[3] 李震東, 趙建福, 魯仰輝, 等. 池沸騰現(xiàn)象中熱毛細(xì)對(duì)流的成因[J]. 空間科學(xué)學(xué)報(bào), 2008, 28(1): 38-43.

LI Zhendong, Zhao Jianfu, Lu Yanghui, et al. Origin of thermocapillary convection in pool boiling[J]. Chinese Journal of Space Science, 2008, 28(1):38-43.

[4] Tanasawa I. Advances in condensation heat transfer[J]. Advances in Heat Transfer, 1991, 21: 55-139.

[5] 劉秋生, 汪洋, 紀(jì)巖. 蒸發(fā)相變與界面流動(dòng)耦合機(jī)理研究[J]. 工程熱物理學(xué)報(bào), 2010(10): 1751-1754.

Liu Qiusheng, Wang Yang, Ji Yan. Coupling mechanism of evaporation phase-change and interfacial flow[J]. Journal of Engineering Thermophysics, 2010(10): 1751-1754.

[6] 程向華, 厲彥忠, 陳二鋒, 等. 新型運(yùn)載火箭射前預(yù)冷液氧貯箱熱分層的數(shù)值研究[J]. 西安交通大學(xué)學(xué)報(bào), 2008, 42(9): 1132-1136.

Cheng Xianghua, Li Yanzhong, Chen Erfeng, et al. Numerical investigation of thermal stratification in liquid Oxygen tank for new-style launch vehicle during ground precooling[J]. Journal of Xi'an Jiaotong University, 2008, 42(9): 1132-1136.

[7] Zilliac G, Karabeyoglu M A. Modeling of propellant tank pressurization[R]. AIAA, 2005-3549, 2005.

[8] Kartuzova O, Kassemi M, Agui J H, et al. A CFD model for the Multipurpose Hydrogen Test Bed (MHTB) ground-based self-pressurization and pressure control experiments[C]. Atlanta: 11th AIAA/ASME Joint Thermophysics and Heat Transfer Conference, 2014.

[9] Knoll R H, Smolak G R. Weightlessness experiments with liquid hydrogen in aerobee sounding rockets, uniform radiant heat addition-flight 1[R]. NASA TM X-484, 1962.

[10] ANSYS FLUENT 13.0 Theory Guide[M]. Canonsburg, PA, ANSYS. Inc., 2011.

Computational Research on Phase Change Model forCryogenic Propellant in Microgravity

Wang Xi1, Wang Jue2, Rong Yi1, Huang Hui1

(1. Beijing Institute of Astronautical Systems Engineering, Beijing, 100076; 2. China Academy of Launch Vehicle Technology, Beijing, 100076)

Boil off is a characteristic of cryogenic propellant, and the phenomenen of phase change is the major factor of cryogenic propellant storage on orbit. As the phase change model plays a key role in cryogenic fluid Computational Fluid Dynamics(CFD) simulation, the usage of reasonable phase change model will be important. Based on four different phase change theory, the commercial software FLUENT are used and four different phase change models are presented by user defined secondary developing code. A CFD research on liquid hydrogen evaporation in microgravity is present, comparing with NASA sounding rocket experiment. Finally, a reasonable phase change model is commended by comparing the pressurizing rate of the four models, which can be used in cryogenic fluid management simulation.

Cryogenic propellant; Phase change; Microgravity

1004-7182(2018)01-0036-05

10.7654/j.issn.1004-7182.20180107

V511

A

2016-12-14；

2017-02-20

王 夕（1989-），男，博士，工程師，主要研究方向?yàn)檫\(yùn)載火箭總體設(shè)計(jì)