李佳超，梁國柱

(北京航空航天大學宇航學院，北京 100083)

0 引 言

低溫火箭貯箱壁面的漏熱使推進劑汽化，引起推進劑的損失和貯箱內(nèi)壓力上升，對貯箱工作過程和工作安全都有不利影響。因此，必須對低溫火箭貯箱的自增壓過程進行研究。

國內(nèi)外對低溫推進劑貯箱的自增壓過程開展過一系列實驗和理論研究。Aydelott[1]對直徑9英寸的球形液氫貯箱進行地面測試，發(fā)現(xiàn)增壓速率受加熱方式布置的影響，但是受充填率變化的影響很小，且增壓速率隨熱量變化近似線性上升。Aydelott和Spuckler[2]研究不同貯箱尺寸的影響，發(fā)現(xiàn)在相同加熱速率和體積比例情況下具有相似的壓力變化速率，且體積充填率對壓力上升速率呈正相關。文獻[3-4]采用4.95 m3的多層絕熱液氫貯箱進行自增壓實驗，給出壁面溫度和充填率對液氫貯箱內(nèi)壓力變化和溫度分布的影響，并發(fā)現(xiàn)自增壓速率隨熱載荷增加而上升，但體積充填率與壓力上升速率并不是嚴格的正相關。Joseph等[5]利用不同絕熱層厚度貯箱研究絕熱層厚度對壓力上升和溫度分層的影響，結(jié)果表明絕熱層厚度的減少會導致壓力上升速率變快，溫度分層加劇。而國內(nèi)直接采用液氫和液氧開展自增壓實驗尚未見報道。王貴仁[6]用液氮作為工質(zhì)，測定了不同初始充填率下壓力隨時間的變化。喬國發(fā)[7]采用液化天然氣(LNG)作為工質(zhì)，給出了密閉貯罐內(nèi)壓力、日蒸發(fā)氣體量及蒸發(fā)率間關系。以上實驗多注重于工程方面的應用，缺乏對自增壓過程實驗數(shù)據(jù)及規(guī)律的研究，也沒有為數(shù)值仿真提供準確的初邊值條件。自增壓的數(shù)值模擬主要有集總參數(shù)法和計算流體力學(CFD)法，如文獻[8]將貯箱內(nèi)流體分為氣枕區(qū)、液體區(qū)和飽和區(qū)，飽和區(qū)作為傳熱傳質(zhì)的交換區(qū)。文獻[9]在三分區(qū)法的基礎上增加了邊界層區(qū)，用于考慮壁面附近液體的流動。Brasi和Kassemi[10]對氣枕區(qū)采用集總參數(shù)法，液體區(qū)采用CFD方法，通過與液體區(qū)輸運方程耦合求解氣枕區(qū)集總參數(shù)。分區(qū)法的缺陷是不同分區(qū)間的換熱系數(shù)難以準確地確定。為了彌補集總參數(shù)法的不足， Liu等[11]運用CFD方法對液氧貯箱建立二維模型，給出了貯箱內(nèi)溫度分布和壓力的變化，并且通過對比得出二維模型在貯箱增壓過程中具有更高的準確性。

本文采用液氮作為模擬介質(zhì)，設計了可視化液氮貯箱自增壓實驗系統(tǒng)，通過實驗測量貯箱內(nèi)壓力、溫度和液面的變化，分析初始充填率與貯箱內(nèi)壓力的變化關系。基于實驗現(xiàn)象和測定的初邊值條件，選擇合理的物理模型和CFD方法對液氮貯箱自增壓過程進行數(shù)值仿真，深入分析整個流場的細節(jié)與機理，為低溫推進劑貯箱的自增壓過程提供合理的預測模型。

1 實驗系統(tǒng)

為了給仿真提供準確的初邊值條件，并觀測低溫流體流動和液面變化情況，建立了液氮貯箱自增壓實驗系統(tǒng)，如圖1所示。

實驗貯箱采用高純度石英玻璃，貯箱總高478 mm，圓柱段內(nèi)徑290 mm，出口內(nèi)徑98 mm，平均厚度8 mm，總體積23.71 L。貯箱自增壓實驗進行前需要對整個系統(tǒng)的密封性進行測試，在0.2 MPa(絕壓，下同)下經(jīng)過多次壓力實驗得到系統(tǒng)的平均漏氣率為5 Pa/s。

2 實驗結(jié)果與分析

2.1 自增壓壓力變化規(guī)律

圖2給出了貯箱內(nèi)壓力、溫度隨時間和體積充填率的變化關系。根據(jù)實驗測得的壓力數(shù)據(jù)認為貯箱內(nèi)液氮自增壓過程的壓力變化可劃分為具有典型意義的初始段、過渡段和穩(wěn)定段，在圖中分別以數(shù)據(jù)1,2,3表示。從圖2可以看出，增壓過程貯箱內(nèi)液體從過飽和狀態(tài)向過冷狀態(tài)進行轉(zhuǎn)變，初始充填率越大，液體的過冷度越大。液體主體區(qū)(液體表層區(qū)和壁面附近以外的液體區(qū)域)溫度僅在初始段短時間內(nèi)高于氣枕區(qū)壓力對應的飽和溫度，在初始段其余時間內(nèi)及整個過渡段、穩(wěn)定段，液體主體區(qū)溫度都低于對應氣枕區(qū)壓力下的飽和溫度。初始段是一個短暫的壓力迅速上升過程，原因是：自增壓開始時，液體的溫度略高于氣枕區(qū)壓力下對應的飽和溫度，液體處于劇烈的沸騰狀態(tài)，液體汽化的質(zhì)量迅速增加，且汽化的蒸汽溫度較低，貯箱內(nèi)壁面對低溫蒸汽有劇烈的加熱作用。過渡段是自增壓中壓力緩慢上升的過程，實驗中觀測到過渡段液體由劇烈沸騰過渡到貯箱內(nèi)壁面局部沸騰，內(nèi)部液體運動速度緩慢，壁面處有小氣泡向界面處運動，氣液間的對流運動受到抑制，使氣液界面變?yōu)榉€(wěn)定的水平面。原因是：初始段壓力的迅速上升使得液體汽化量迅速下降，甚至出現(xiàn)氣體的冷凝，貯箱內(nèi)的液體處于過冷狀態(tài)，外部的熱量大部分用于液體的升溫。穩(wěn)定段在自增壓中持續(xù)時間最長，其壓力上升速率介于初始段和過渡段的壓力上升速率之間，液體的汽化位置集中于壁面附近的液體區(qū)域和液體表層區(qū)(氣液界面附近區(qū)域)，氣液界面是穩(wěn)定的水平面。原因是：隨著時間的增加，液體主體區(qū)的過冷度增大，氣枕區(qū)溫度升高，與貯箱內(nèi)壁面的換熱作用減弱，壓力上升主要來自液體表層區(qū)和貯箱壁面附近液體的汽化。壓力上升速率趨于穩(wěn)定的原因是貯箱內(nèi)液體的汽化速率趨于穩(wěn)定。綜合圖2可知:初始段壓力上升的主要因素是液體的劇烈沸騰和貯箱壁面對氣枕區(qū)的加熱作用，過渡段壓力上升的主要因素是貯箱壁面對氣枕區(qū)的加熱作用，而穩(wěn)定段壓力上升的主要因素是液體的汽化。

圖3給出了貯箱內(nèi)氣枕區(qū)壓力與體積充填率的關系，從圖3(a)可以看出，在低體積充填率下，體積充填率越大，同一時刻貯箱內(nèi)的壓力越大。初始段和穩(wěn)定段的壓力上升速率隨體積充填率增大而增大。而過渡段在低體積充填率下不明顯。從圖3(b)可以看出，在高體積充填率下，過渡段的壓力上升速率低，持續(xù)時間長。穩(wěn)定段的壓力上升速率隨體積充填率增大而增大。從圖3(c)可以看出，穩(wěn)定段壓力上升速率在體積充填率為53.5%時明顯高于體積充填率為45.2%時。綜合圖3可知：穩(wěn)定段的壓力上升速率隨體積充填率的增大而增大；過渡段在低體積充填率下作用不明顯，高體積充填率下持續(xù)時間長；初始段的壓力上升速率與體積充填率間存在一個極大值，使初始段能達到最大的壓力值。實驗中觀測到液體的汽化主要在貯箱內(nèi)壁面和液體表層區(qū)。初始段氣枕區(qū)溫度低，壁面的加熱作用大，液體的汽化量大，體積充填率的增加使固體壁面與氣枕區(qū)的換熱量減少，液體區(qū)與固體壁面換熱量增加。在過渡段，部分氣體冷凝，壓力上升主要取決于固體壁面對氣枕區(qū)的加熱，液體體積充填率越大，液體表層區(qū)需要吸收更多熱量達到飽和狀態(tài)，則持續(xù)時間越長。穩(wěn)定段氣體溫度持續(xù)升高導致氣枕區(qū)與固體壁面的換熱呈下降趨勢，液體的汽化對壓力上升起主導作用，而體積充填率越大，液體汽化量越大。

2.2 自增壓固體區(qū)溫度變化規(guī)律

由于自增壓中液面幾乎不變，圖4只給出初始的氣液界面位置，用虛線表示；細實線表示貯箱不同區(qū)域的劃分。從圖4(a)可以看出，氣液界面以下測點的溫度隨時間變化緩慢，氣液界面以上測點溫度隨時間有較大的變化(上部不銹鋼封頭處測點除外，因其導熱率遠高于石英玻璃)，貯箱外壁面溫度在軸向方向存在溫度梯度，距出口位置越近溫度越高。從圖4(b)可以看出，50 s之前氣液界面的貯箱柱狀段外壁面測點溫度偏低，50 s后液體區(qū)的貯箱外壁面溫度在軸向存在溫度梯度。從圖4(c)可以看出，液體區(qū)的貯箱外壁面溫度在靠近氣液界面偏低，沒有明顯的溫度梯度分布。綜合圖4得出，自增壓過程中貯箱氣枕區(qū)固體外壁面溫度高于液體區(qū)固體外壁面溫度，氣枕區(qū)固體外壁面溫度隨時間變化速度快，液體區(qū)固體壁面溫度隨時間變化速度慢。液體的體積充填率越高，貯箱外壁面的溫度越低。原因是：氣體密度小，在比熱容差距不大的情況下，氣枕區(qū)的溫度上升速率快，使得外界傳遞的熱量大部分用于固體區(qū)域的升溫，導致氣枕區(qū)貯箱外壁面溫度高。高體積充填率氣液界面附近溫度偏低是因為氣液界面附近液體的汽化速度快，使得液體與固體區(qū)熱量交換值大，外界用于固體區(qū)溫度上升的熱量變小。貯箱外壁面溫度分布說明了貯箱與外界熱交換值與貯箱的形狀和液體充填率有關。

2.3 自增壓氣枕區(qū)和液體區(qū)溫度變化規(guī)律

圖5給出了流體溫度隨時間的變化關系。從圖5可以看出：貯箱內(nèi)氣枕區(qū)存在較大的溫度梯度，氣枕區(qū)溫度高于液體區(qū)溫度，使得氣枕區(qū)向液體表層區(qū)傳熱，并且體積充填率越大，氣枕區(qū)的溫度越小。液體表層區(qū)溫度高，液體主體區(qū)域各溫度測點相差不大，且體積充填率越大，液體區(qū)溫度越趨于一致。貯箱內(nèi)氣枕區(qū)和液體區(qū)的溫度隨時間的增加而增大，并且氣枕區(qū)溫度增速高于液體區(qū)溫度的增速。液體區(qū)溫度分層原因是：在固體壁面和氣枕區(qū)的加熱作用下，壁面附近液體和液體表層區(qū)溫度升高，壁面附近溫度較高的液體在浮力驅(qū)動下向上運動，使得液體溫度趨于一致，而液體表層區(qū)溫度升高抑制了液體表層區(qū)和液體主體區(qū)的對流運動，隨著時間的推移，浮力驅(qū)動液體運動作用減弱，液體表層溫度上升抑制液體的對流運動作用增強，液體逐漸趨于靜止，這與實驗中加壓后液體幾乎靜止現(xiàn)象吻合。而氣枕區(qū)溫度隨體積充填率增大而減小的原因是：液體區(qū)與固體壁面接觸面積越大，使得壁面附近汽化生成的氣體越多，液體汽化生成氣體的溫度比氣枕區(qū)溫度低，二者的混合導致氣枕區(qū)溫度偏低。

3 數(shù)值仿真與實驗結(jié)果比較

3.1 數(shù)學模型

實驗中觀測到液氮的常壓停放是一個液體劇烈沸騰過程，而液氮的自增壓是一個液體緩慢運動的過程。由于實驗中對貯箱內(nèi)流體的溫度分布和體積分數(shù)分布缺乏精確的測量，因此通過常壓停放數(shù)值仿真給出自增壓開始時貯箱內(nèi)氣液兩相的溫度和體積分數(shù)分布。

3.1.1控制方程

自增壓的數(shù)值模擬多采用文獻[11-12]所用的流體體積(Volume of fluid,VOF)模型。VOF 模型適合于求解兩種及多種互不相融、分界面相對明顯的多相流問題。對于本文研究自增壓過程，兩相之間基本沒有相互融合、摻混，運動相對獨立，因此選用VOF模型。流體控制方程中的連續(xù)性方程、動量方程和能量方程見文獻[13]。連續(xù)性方程和能量方程考慮液氮和氮氣間傳質(zhì)過程引起質(zhì)量和能量的變化，動量方程中不考慮傳質(zhì)引起的變化。

3.1.2相變模型

相變模型主要有相變系數(shù)模型、假定界面溫度模型、能量平衡模型和溫度恢復模型，后三種模型計算時容易出現(xiàn)不穩(wěn)定，并且需要初始化氣液界面。Knudsen[14]根據(jù)氣體動力學理論推導的相變系數(shù)模型為：

(1)

(2)

式中：下標l,v,vl和lv分別代表液體、氣體、汽化和冷凝；β為控制相變強弱的時間松弛系數(shù)， 單位s-1；α為體積分分數(shù)；ρ為密度， 單位kg/m3；T為溫度， 單位K。J為氣液界面相變質(zhì)量， 單位kg/(m3·s)；傳質(zhì)過程引起的能量交換為式(3)，其中,Q為相變引起的換熱量， 單位J/(m3·s)；ΔH為氣液的相變潛熱， 單位J/kg。

Q=J·ΔH

(3)

相變系數(shù)模型關鍵是確定時間松弛系數(shù)，過大的時間松弛系數(shù)使計算難以收斂，過小的時間松弛系數(shù)會使氣液界面溫度與飽和溫度存在較大偏差[15]。Lee[16]采用固定的時間松弛系數(shù)0.1 s-1取得很好的計算結(jié)果，因此，本文時間松弛系數(shù)取0.1 s-1。

3.1.3固體壁面導熱控制方程

貯箱固體壁面與內(nèi)外流體間存在熱交換作用。固體壁面的導熱方程如下：

(4)

式中：下標w為固體，h為顯焓，單位J/kg；k為導熱系數(shù)，單位W/(m·K)。導熱系數(shù)隨溫度的變化較大，具體值由文獻[17]數(shù)據(jù)采用一階插值計算得出。

3.2 計算方法

貯箱的仿真模型通過Fluent建立，分別對氣枕區(qū)、液體區(qū)和固體區(qū)劃分結(jié)構(gòu)化網(wǎng)格。網(wǎng)格在靠近貯箱內(nèi)壁面和氣液界面的區(qū)域進行適當加密,通過網(wǎng)格數(shù)目無關性驗證計算，最終采用網(wǎng)格數(shù)為37359的二維軸對稱網(wǎng)格，如圖6所示。

貯箱常壓停放過程的初始條件設為氣液兩相飽和溫度77.35 K，固體區(qū)域初始溫度77.35 K，環(huán)境溫度為298 K，貯箱壁面溫度采用圖4(b)中t=0 s測量值。采用恒定壓力出口條件，其值101325 Pa，液氮的體積充填率43.7%。計算時采用VOF兩相流模型，湍流模型采用標準k-ε模型，壁面采用增強壁面函數(shù)，相變模型通過UDF(User-define function)以源項的形式加載到連續(xù)性方程和能量方程。選擇基于壓力的求解器，壓力插值方式選擇Body-force- weight方案，速度壓力耦合方式選擇PISO算法格式，體積分數(shù)插值方式選擇Quick格式，其余的選擇二階迎風格式。氮氣采用理想氣體模型，液氮的密度采用泊松模型。其它物性參數(shù)考慮隨溫度的變化，通過查找文獻[18]采用一階插值計算獲得。時間步長在0.001～0.005 s時，計算結(jié)果基本一致，實際計算時間步長取0.0025 s。通過常壓停放階段的仿真，獲得的內(nèi)外壁面的換熱量如圖7所示，可以看出，165 s后內(nèi)外壁面換熱量近似相等，達到準穩(wěn)定狀態(tài)，因此，選擇此時的狀態(tài)作為自增壓過程的初始條件。此時體積充填率為42.5%，壓力出口調(diào)整為恒定溫度壁面邊界條件。計算方法和常壓停放階段相同，時間步長取0.001 s。

3.3 仿真與實驗結(jié)果比較

圖8給出自增壓實驗與仿真壓力變化的比較。可以看出，仿真中壓力上升曲線與實驗中壓力曲線在變化趨勢相同。仿真中30 s左右后進入穩(wěn)定段，壓力在0.16 MPa左右；而實驗中只用了20 s左右，壓力在0.14 MPa左右。實驗中20 s后的壓力平均上升速率為390 Pa/s，仿真中30 s后壓力平均上升速率為620 Pa/s；穩(wěn)定段的壓力上升速率仿真與實驗的比值在1.58左右。仿真中壓力值偏大的原因：一是液體在汽化過程中需要一定的過熱度，而仿真模型認為液體達到飽和溫度即汽化，造成液體汽化生成的氣體量偏大;二是實驗中未能測量低溫條件下石英玻璃的導熱系數(shù)具體值，而是通過文獻[18]給出的值采用一階插值法獲得，與實際值存在一定的偏差。

圖9給出了仿真得到的貯箱內(nèi)壁面與氣液間的熱交換量?？梢钥闯?，貯箱內(nèi)壁面與液體的換熱量先下降再上升后趨于穩(wěn)定，氣體與壁面的換熱量先下降后趨于穩(wěn)定，穩(wěn)定段液體的換熱量是氣體換熱量的7倍左右。由此可以得出在已有文獻[4, 10, 19-20]的仿真中采用平均熱流密度的邊界條件會使壓力上升速率低于實驗值，與實際物理問題存在較大的偏差。

液氮貯箱自增壓仿真過程的氣液兩相分布及相變位置如圖10所示。圖中黑色區(qū)為液氮汽化區(qū)域，白色區(qū)為氮氣冷凝區(qū)，灰色區(qū)為沒有發(fā)生相變區(qū)。從圖10可以看出：t=0 s時刻，貯箱內(nèi)的液體絕大部分處于沸騰狀態(tài)，液氮汽化為氮氣，只有極小部分氮氣冷凝。t= 150 s時，液氮的汽化主要發(fā)生在貯箱內(nèi)壁面附近和液體表層區(qū)，液體區(qū)有氣體的冷凝，且液體區(qū)內(nèi)冷凝區(qū)域隨時間的推移而擴大。冷凝區(qū)擴大的原因是：飽和溫度上升速度高于液體區(qū)溫度的上升速度，導致液體區(qū)的過冷度隨時間推移而增大。

液氮貯箱自增壓仿真的溫度及速度矢量分布如圖11和圖12所示。從圖11可以看出：氣枕區(qū)沿軸向存在溫度分層，氣體越靠近貯箱的出口溫度越高，而徑向方向同一高度處，靠近壁面區(qū)域的氣體溫度高，氣枕區(qū)的溫度隨時間推移而增加。而液體區(qū)溫度分布隨時間變化顯著不同，t=0 s時，液體區(qū)沿徑向方向表現(xiàn)出溫度分層，從液體內(nèi)部核心區(qū)域向外依次增加，靠近壁面和氣枕區(qū)的溫度高，t=150 s時，液體區(qū)在徑向方向溫度分層消失，液體軸向溫度分層增強，液體區(qū)溫度和軸向溫度梯度也增加，且液體在靠近氣枕區(qū)溫度高。從圖12可以看出：t=0 s時，壁面附近和軸線區(qū)域附近液體向上運動，液體在氣液界面附近形成了順時針的漩渦，氣枕區(qū)在壁面附近沿壁面向上運動，氣液界面附近氣體同液體的運動趨勢相同。t=150 s時，壁面附近流體的運動同t=0 s時相同，液體區(qū)向上運動得到抑制，液體向遠離氣液界面處運動。原因是：液體在氣枕的加熱作用下，上部溫度升高抑制了液體的對流運動，液體的整體溫度在氣枕和壁面加熱作用下逐漸上升。氣枕區(qū)在氣液界面附近向液體區(qū)傳遞能量，導致溫度偏低，并且在壁面加熱作用下，氣枕區(qū)整體溫度升高，由于氣枕區(qū)的熱容量小，溫度上升速度遠高于液體區(qū)。

4 結(jié) 論

本文對液氮貯箱自增壓過程進行了實驗和仿真，得出如下結(jié)果：

1)貯箱內(nèi)液氮自增壓過程初始段的壓力上升速率最快，過渡段的壓力上升速率最小，穩(wěn)定段的壓力上升速率近似保持穩(wěn)定。

2)貯箱內(nèi)自增壓過程壓力變化曲線中，穩(wěn)定段的壓力上升速率隨液體體積充填率的增大而增大，初始段的壓力上升速率取決于液體的過熱度，而過渡段持續(xù)時間與體積充填率有關。

3)貯箱內(nèi)液氮自增壓過程固體區(qū)域外壁面溫度在低體積充填率下沿軸向線性分布，在高體積充填率下氣液界面附近的測點溫度顯著偏低。

4)貯箱內(nèi)液氮自增壓是一個氣液界面呈現(xiàn)為靜止的水平面，液體對流受到抑制，氣枕區(qū)和液體區(qū)存在溫度分層的過程。

5)CFD仿真采用相變系數(shù)模型模擬貯箱從常壓停放到自增壓整個過程，能夠比較準確地預測貯箱內(nèi)壓力和溫度變化規(guī)律、液面行為、相變位置等。

[1] Aydelott J C. Normal gravity self-pressurization of 9-inch (23 cm) diameter spherical liquid hydrogen tankage [R].Cleveland, USA：Lewis Research Center, 1967.

[2] Aydelott J C, Spuckler C M. Effect of size on normal-gravity self-pressurization of spherical liquid hydrogen tankage[R]. Cleveland, USA: Lewis Research Center, 1969.

[3] Hasan M M, Lin C S, Vandresar N T. Self-pressurization of a flight weight liquid hydrogen storage tank subjected to low heat flux [C].ASME/AICHE National Heat Transfer Conference, Minneapolis, USA, 1991.

[4] Stochl R J, Knoll R H. Thermal performance of a liquid hydrogen tank multilayer insulation system at warm boundary temperatures of 630, 530, and 152 R [C]. The 27th Joint Propulsion Conference, Sacramento,USA, 1991.

[5] Joseph J, Agrawal G, Agarwal D K, et al. Effect of insulation thickness on pressure evolution and thermal stratification in a cryogenic tank [J].Applied Thermal Engineering, 2017, 111: 1629-1639.

[6] 王貴仁. 低溫容器內(nèi)壓力變化規(guī)律的研究[D].蘭州：蘭州理工大學, 2008.[Wang Gui-ren. The research on variation law of pressure in cryogenic vessel [D].Lanzhou: Lanzhou University of Technology, 2008.]

[7] 喬國發(fā). 影響LNG儲存容器蒸發(fā)率因素的研究[D].東營：中國石油大學, 2007.[Qiao Guo-fa. The study on the influential factors of evaporation rate of the liquefied natural gas tank [D]. Dongying: China University of Petroleum, 2007.]

[8] Zilliac G, Karabeyoglu M A. Modeling of propellant tank pressurization[C].The 41st Joint Propulsion Conference & Exhibit,Tucson,USA,2005.

[9] 劉展, 厲彥忠, 王磊. 低溫推進劑熱分層研究[J]. 宇航學報, 2015, 36(6): 613-623. [Liu Zhan, Li Yan-zhong, Wang Lei. Research on cryogenic propellant thermal stratification [J]. Journal of Astronautics, 2015, 36(6): 613-623.]

[10] Barsi S, Kassemi M. Numerical and experimental comparisons of the self-pressurization behavior of an LH2 tank in normal gravity [J]. Cryogenics, 2008, 48(3): 122-129.

[11] Liu Z, Li Y, Jin Y. Pressurization performance and temperature stratification in cryogenic final stage propellant tank [J]. Applied Thermal Engineering, 2016, 106: 211-220.

[12] Chen L, Liang G. Simulation research of vaporization and pressure variation in a cryogenic propellant tank at the launch site [J]. Microgravity Science and Technology, 2013, 25(4): 203-211.

[13] Fluent A. Ansys fluent theory guide[EB/OL].2010[2017].https://support.ausys.com/.

[14] Knudsen M. Kinetic theory of gases: some modern aspects [M].London: Methuen, 1950: 103-160.

[15] De Schepper S C K, Heynderickx G J, Marin G B. Modeling the evaporation of a hydrocarbon feedstock in the convection section of a steam cracker [J]. Computers & Chemical Engineering, 2009, 33(1): 122-132.

[16] Lee W H. Pressure iteration scheme for two-phase flow modeling [R]. Los Alamos,USA:Los Alamos National Laboratory,1979.

[17] 陳國邦. 低溫工程材料[M]. 杭州: 浙江大學出版社, 1998: 195-201.

[18] 陳國邦. 低溫流體熱物理性質(zhì)[M]. 北京：國防工業(yè)出版社, 2006: 190-220.

[19] Panzarella C H, Kassemi M. On the validity of purely thermodynamic descriptions of two-phase cryogenic fluid storage [J]. Journal of Fluid Mechanics, 2003, 484: 41-68.

[20] Mattick S J, Lee C P, Hosangadi A, et al. Progress in modeling pressurization in propellant tanks [C]. Joint Propulsion Conference & Exhibit, Nashville, USA, 2010.