零重力條件下低溫射流抑制大尺寸液氫儲(chǔ)罐熱分層的數(shù)值研究1)

郭 斌 趙建福 李 凱,3) 胡文瑞

?(中國(guó)科學(xué)院大學(xué)工程科學(xué)學(xué)院，北京 100049)

?(中國(guó)科學(xué)院力學(xué)研究所微重力重點(diǎn)實(shí)驗(yàn)室，北京 100190)

??(中國(guó)科學(xué)院力學(xué)研究所高溫氣體動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京 100190)

引言

以液氫為代表的低溫推進(jìn)劑在航天領(lǐng)域有著廣泛的應(yīng)用，但是因?yàn)槠浞悬c(diǎn)低、易蒸發(fā)的特性，在地面停放、發(fā)射過程以及在軌運(yùn)行等各個(gè)階段都極易發(fā)生汽化，引起儲(chǔ)罐內(nèi)部熱分層和自增壓，進(jìn)而威脅飛行器系統(tǒng)結(jié)構(gòu)安全.在軌運(yùn)行期間，推進(jìn)劑儲(chǔ)罐會(huì)受到各種空間輻射的影響，比如太陽輻射、地球反照熱外流、地球紅外輻射熱流以及空間黑背景輻射熱流等等.這些熱源的存在都會(huì)使得推進(jìn)劑儲(chǔ)罐因受熱不均勻在局部發(fā)生汽化現(xiàn)象.傳統(tǒng)的排氣方法能有效降低儲(chǔ)罐內(nèi)部壓強(qiáng)，但一定程度上造成了推進(jìn)劑的浪費(fèi)，使發(fā)射成本增加，同時(shí)制約飛行器運(yùn)載能力，不利于長(zhǎng)期在軌飛行及深空探測(cè)活動(dòng).20 世紀(jì)末，美國(guó)NASA 的研究者提出了零蒸發(fā)存儲(chǔ)(zero boil-off,ZBO)的概念.ZBO 技術(shù)是通過將被動(dòng)絕熱技術(shù)與主動(dòng)制冷技術(shù)有機(jī)結(jié)合起來，實(shí)現(xiàn)低溫推進(jìn)劑無損儲(chǔ)存的技術(shù)，其中被動(dòng)絕熱技術(shù)主要是使用絕熱材料減少外部熱量的滲入.Hastings 等[1]發(fā)現(xiàn)對(duì)于需要長(zhǎng)期在軌運(yùn)行的儲(chǔ)罐而言，單純使用被動(dòng)絕熱技術(shù)會(huì)增加儲(chǔ)罐質(zhì)量，進(jìn)而降低航天器的靈活性，因此在使用絕熱材料的基礎(chǔ)上嘗試主動(dòng)制冷技術(shù)的運(yùn)用對(duì)于長(zhǎng)期航天探索工程具有重要意義.主動(dòng)制冷技術(shù)主要有兩種，一種是通過熱傳導(dǎo)裝置將罐內(nèi)熱量傳導(dǎo)到罐外，另一種則是通過機(jī)械攪拌或者低溫射流的方式加強(qiáng)儲(chǔ)罐內(nèi)部流體對(duì)流.熱傳導(dǎo)裝置作用的范圍相對(duì)較小，機(jī)械攪拌的實(shí)現(xiàn)工藝更加復(fù)雜，因此通過施加儲(chǔ)罐內(nèi)部低溫射流的方式進(jìn)行壓強(qiáng)控制成為國(guó)內(nèi)外相關(guān)領(lǐng)域?qū)W者的研究熱點(diǎn).

通過低溫射流實(shí)現(xiàn)熱分層和自增壓控制的研究可以追溯到20 世紀(jì)70 年代，Poth 等[2]對(duì)比分析了各種熱分層消除裝置的特點(diǎn)，認(rèn)為軸向低溫射流具有很高的流體混合性能并且輕便靈活，增壓控制易于實(shí)現(xiàn)；Merte 等[3]研究了軸對(duì)稱的圓柱形儲(chǔ)罐中單組分工質(zhì)氣相和液相之間的傳熱和傳質(zhì)過程，其中氣相采用理想氣體假設(shè)，將控制方程轉(zhuǎn)換為有限差分形式，進(jìn)行數(shù)值計(jì)算得到了多種計(jì)算工況下的壓力演化曲線.Audelott[4]在Lewis Research Center 5-10 s零重力裝置中實(shí)驗(yàn)研究了向直徑為10 cm的圓柱形儲(chǔ)罐中軸向噴射乙醇時(shí)產(chǎn)生的液體流動(dòng)模式，發(fā)現(xiàn)儲(chǔ)罐內(nèi)部流場(chǎng)分布與儲(chǔ)罐幾何形狀、低溫射流速度、儲(chǔ)罐填充比以及低溫射流位置有關(guān)；Lin 等[5]先后實(shí)驗(yàn)研究了體積為0.144 m3的小型儲(chǔ)罐和體積為4.89 m3大型儲(chǔ)罐使用低溫射流混合裝置實(shí)現(xiàn)增壓控制的效果，實(shí)驗(yàn)表明低溫射流可以有效控制罐內(nèi)增壓，但是小罐模型得出的低溫射流時(shí)間和增壓速率之間的關(guān)系無法很好預(yù)測(cè)大罐模型中兩者之間的相對(duì)關(guān)系；Panzarella 等[6]研究了常重力條件下小型儲(chǔ)罐內(nèi)部氣液兩相的自增壓情況，證實(shí)了增強(qiáng)儲(chǔ)罐內(nèi)部對(duì)流可以抑制熱分層現(xiàn)象進(jìn)而實(shí)現(xiàn)增壓控制；Mukka 等[7]研究了常重力條件下低溫射流方式對(duì)消除熱分層效果的影響，認(rèn)為不同低溫射流條件導(dǎo)致的儲(chǔ)罐內(nèi)部流場(chǎng)的差異對(duì)儲(chǔ)罐內(nèi)溫度場(chǎng)和速度場(chǎng)的分布有著重要影響.在認(rèn)識(shí)到低溫射流可以有效消除熱分層現(xiàn)象后，科研工作者建立了很多罐體系統(tǒng)研究如何優(yōu)化低溫射流消除熱分層的效果.Ho 等[8]建立了一個(gè)泵--熱管模型，通過泵裝置促進(jìn)儲(chǔ)罐內(nèi)部流體循環(huán)，不斷把高溫流體噴向熱管，再由熱管將熱量導(dǎo)出.通過對(duì)該模型的三維數(shù)值模擬發(fā)現(xiàn)通過增大低溫射流速度可以有效消除儲(chǔ)罐內(nèi)部的高溫區(qū)域；隨后Ho 等[9]利用數(shù)值模擬方法研究了軸向圓盤多孔低溫射流噴嘴結(jié)構(gòu)，發(fā)現(xiàn)在保持入射速度不變的情況下，加大孔口直徑也可以顯著提升熱分層消除效果，此外射流噴嘴位置對(duì)熱分層消除效果有著顯著影響,并對(duì)圓盤射流噴嘴系統(tǒng)進(jìn)行了參量分析，研究了入口管直徑、射流噴嘴位置和射流噴嘴直徑對(duì)系統(tǒng)制冷效果的影響.Ho 等[10]又對(duì)之前建立的泵--熱管模型進(jìn)行了深入研究，分析了噴射間隔和噴管長(zhǎng)度等參量對(duì)消除熱分層效果的影響.Belmedani 等[11]開展了不同熱通量的液氮儲(chǔ)罐熱分層的實(shí)驗(yàn)研究，熱通量越大，自由面蒸發(fā)速率越大，自由面和液相區(qū)溫差越明顯，較好地揭示了蒸發(fā)機(jī)理，通過研究低溫儲(chǔ)罐中液體的速度和溫度特性，認(rèn)為強(qiáng)熱流誘發(fā)熱分層現(xiàn)象.Zilliac 等[12]提出了一種平衡熱力學(xué)模型來預(yù)測(cè)儲(chǔ)罐增壓過程，該模型適用于具有高揮發(fā)性的低溫推進(jìn)劑，通過與實(shí)測(cè)結(jié)果進(jìn)行對(duì)比，證實(shí)了該模型的準(zhǔn)確性.Grayson 等[13]模擬了常重力情況下通過外部加熱實(shí)現(xiàn)增壓和通過熱力學(xué)排氣實(shí)現(xiàn)減壓的儲(chǔ)罐數(shù)值計(jì)算模型，假設(shè)液相不可壓縮，其密度僅是溫度的函數(shù).該模型的壓力和溫度預(yù)測(cè)結(jié)果與測(cè)試中的傳感器測(cè)量結(jié)果相比較有很高的一致性.Barsi 等[14]采取了與原有的集總熱力學(xué)模型不同的氣液兩相模型，在該模型中，液相和汽相的守恒方程都得到了求解.研究中為了簡(jiǎn)化計(jì)算將兩相都視為不可壓縮的，并使用積分方法來求解交界面的質(zhì)量交換.通過將儲(chǔ)罐壓力上升的數(shù)值計(jì)算與以往計(jì)算模型的結(jié)果進(jìn)行比較，證實(shí)了該模型的準(zhǔn)確性.Lopez 等[15]提出了一種新的計(jì)算流體動(dòng)力學(xué)模型，模擬了重力作用下橢球形液氫儲(chǔ)罐在外部加熱情形下的壓力控制過程.壓力控制的實(shí)現(xiàn)由位于儲(chǔ)罐內(nèi)部的軸向射流熱動(dòng)力排氣系統(tǒng)(TVS)提供，該系統(tǒng)向儲(chǔ)罐內(nèi)注入冷流體，使液體發(fā)生混合以降低儲(chǔ)罐內(nèi)部壓力.該軸對(duì)稱模型的計(jì)算采用商業(yè)軟件FLOW-3D 進(jìn)行計(jì)算，其中定量模型驗(yàn)證采用1999 年在馬歇爾航天飛行中心進(jìn)行的工程檢驗(yàn)測(cè)試結(jié)果，計(jì)算結(jié)果表明模型預(yù)測(cè)的自增壓速率和流體溫度與試驗(yàn)數(shù)據(jù)吻合良好.這項(xiàng)研究提升了當(dāng)時(shí)用于實(shí)現(xiàn)低溫壓力控制的CFD 建模能力，并為開發(fā)基于CFD 的空間硬件設(shè)計(jì)提供了參考.Barsi 等[16]提出了一種兩相CFD 模型，該模型描述了在常重力條件下，部分填充的LH2 儲(chǔ)罐的自增壓行為，并利用已有的不同填充水平下的實(shí)驗(yàn)數(shù)據(jù)，對(duì)模型的預(yù)測(cè)能力進(jìn)行了評(píng)估.評(píng)估結(jié)果表明該模型的預(yù)測(cè)結(jié)果與實(shí)驗(yàn)測(cè)得的壓力曲線吻合較好.Kumar 等[17]計(jì)算了不同縱橫比的大型液氫儲(chǔ)罐中蒸發(fā)對(duì)儲(chǔ)罐熱分層的影響.他們采用均勻兩相模型，分別求解了氣相和液相各自的守恒方程，氣液界面的蒸發(fā)是通過傳質(zhì)源項(xiàng)來實(shí)現(xiàn)的，隨著儲(chǔ)罐的縱橫比增加，分層的程度會(huì)逐漸增加.Oliveira 等[18]建立了一個(gè)儲(chǔ)罐熱分層模型，其中包括了平臺(tái)的調(diào)節(jié)旋轉(zhuǎn).該模型可用于評(píng)估軸向加速度、自旋速率、熱通量以及儲(chǔ)罐幾何形狀等因素對(duì)推進(jìn)劑儲(chǔ)罐內(nèi)部熱分層現(xiàn)象的影響.Li 等[19]通過實(shí)驗(yàn)和數(shù)值方法研究了多層絕緣低溫儲(chǔ)罐中真空損耗引起的液氮熱分層的瞬態(tài)過程，他們認(rèn)為真空損耗可以迅速導(dǎo)致熱分層.在實(shí)驗(yàn)中，利用熱電偶得到了實(shí)驗(yàn)儲(chǔ)罐中液體溫度的分布和演變.通過建立二維模型進(jìn)行數(shù)值模擬計(jì)算，模擬了儲(chǔ)罐內(nèi)部熱分層的形成和減弱過程以及液體溫度場(chǎng)分布，并將數(shù)值模擬結(jié)果與實(shí)驗(yàn)結(jié)果進(jìn)行了比較.研究還表明，罐內(nèi)的兩相流動(dòng)對(duì)熱分層起著重要作用.Wang 等[20]通過數(shù)值模擬系統(tǒng)對(duì)比了微重力條件下由熱管和噴管組成的ZBO 系統(tǒng)中，噴管數(shù)量、噴管出流方向和熱管蒸發(fā)器熱傳遞效率對(duì)增壓控制性能的影響.Liu 等[21]通過在Rahman 模型的基礎(chǔ)上增加了向下的導(dǎo)管，使儲(chǔ)罐底部也發(fā)生強(qiáng)制對(duì)流，提升ZBO 性能.Liu 等[22]進(jìn)一步通過正交試驗(yàn)設(shè)計(jì)，研究了導(dǎo)管出口到罐底的距離、噴管的半徑、射流噴嘴深度、導(dǎo)管半徑和環(huán)形射流噴嘴直徑等參量對(duì)儲(chǔ)罐設(shè)計(jì)的影響.Roh 等[23]使用商業(yè)軟件FLUENT 對(duì)壓縮液化天然氣儲(chǔ)罐中的瞬態(tài)自然對(duì)流進(jìn)行了數(shù)值計(jì)算.計(jì)算結(jié)果表明，蒸發(fā)氣體的生成在很大程度上取決于儲(chǔ)罐內(nèi)部垂直溫度的分布，而儲(chǔ)罐內(nèi)部垂直方向的溫度分布受儲(chǔ)罐自增壓過程的影響.同時(shí)作者對(duì)儲(chǔ)罐壓力、罐體尺寸和增壓過程對(duì)蒸發(fā)氣體生成的影響進(jìn)行了量化研究.Wang 等[24]采用計(jì)算流體動(dòng)力學(xué)模型，對(duì)液氫儲(chǔ)罐的加壓排放過程進(jìn)行了數(shù)值模擬，將壁面區(qū)域和流體區(qū)域同時(shí)考慮為計(jì)算域，采用低雷諾數(shù)k-ε 模型來處理流體和壁面的熱交換效應(yīng).該模型還考慮了氣液相變，并將數(shù)值計(jì)算結(jié)果與已有實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比，認(rèn)為該CFD 模型在增壓計(jì)算過程中具有良好的適應(yīng)性.通過該模型可以獲得增壓所需氣體量、儲(chǔ)罐內(nèi)部壓力變化過程、儲(chǔ)罐內(nèi)部溫度分布等詳細(xì)特征.作者還分析了相變效應(yīng)和儲(chǔ)罐結(jié)構(gòu)對(duì)增壓性能的影響，計(jì)算結(jié)果表明，氣液相變對(duì)增壓行為影響不大.隨后Wang 等[25]又建立了一種計(jì)算流體力學(xué)模型，該模型可同時(shí)考慮罐內(nèi)的熱交換和外部空氣動(dòng)力加熱，并對(duì)低溫儲(chǔ)罐在排放過程中的瞬態(tài)熱狀態(tài)和增壓性能進(jìn)行了研究.該模型中計(jì)算域不僅包括了流體和罐壁區(qū)域，還包括了泡沫絕熱層區(qū)域.作者將該模型的計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行比較，兩者具有很高的一致性.然后利用該模型對(duì)儲(chǔ)罐的增壓排放過程進(jìn)行了預(yù)測(cè)，得到了其熱力學(xué)行為、增壓行為過程.Daigle 等[26]提出了一種描述液氫低溫儲(chǔ)罐中自然對(duì)流的溫度分層的動(dòng)態(tài)模型，用MATLAB 實(shí)現(xiàn)了一個(gè)通用的低溫儲(chǔ)罐的溫度分層現(xiàn)象模擬程序，該模型可以模擬常重力以及重力增加和減少條件下儲(chǔ)罐內(nèi)部溫度分層結(jié)果.Fu 等[27]采用數(shù)值計(jì)算的方法，對(duì)部分填充的圓柱肋式液氫儲(chǔ)罐在不同肋間距比下的自增壓過程進(jìn)行了研究.利用商業(yè)軟件FLUENT 進(jìn)行計(jì)算，選取了流體體積法和相變模型，并通過自定義函數(shù)對(duì)模型進(jìn)行適當(dāng)修改，建立了用于預(yù)測(cè)儲(chǔ)罐內(nèi)部流體流動(dòng)和傳熱的求解方案.隨后Fu 等[28]對(duì)低溫儲(chǔ)罐在微重力條件下的蒸發(fā)及其對(duì)蒸汽壓的影響進(jìn)行了數(shù)值研究，研究了表面張力、氣泡接觸角和重力等因素的影響.李佳超等[29]以液氮為研究工質(zhì)利用透明玻璃搭建的低溫儲(chǔ)罐自增壓實(shí)驗(yàn)系統(tǒng)，研究了自增壓過程壓力和溫度的變化規(guī)律及填充比對(duì)壓力和溫度變化過程的影響.2018 年王夕等[30]使用FLUENT 軟件進(jìn)行數(shù)值模擬，對(duì)比研究了4 種相變模型對(duì)微重力環(huán)境中液氫推進(jìn)劑受熱蒸發(fā)過程的影響，并與國(guó)外探空火箭試驗(yàn)數(shù)據(jù)進(jìn)行比較，分析對(duì)比了4 種相變模型各自的優(yōu)劣.李鵬等[31]以液氫、液氧等低溫推進(jìn)劑為研究對(duì)象，開展控制低溫推進(jìn)劑在儲(chǔ)罐內(nèi)部因受熱而蒸發(fā)的現(xiàn)象，通過對(duì)不同條件下的蒸發(fā)量控制效果的對(duì)比分析，為今后航天器推進(jìn)劑儲(chǔ)罐的設(shè)計(jì)提供參考依據(jù).郭志釩等[32]分析了高壓儲(chǔ)氫、低溫液態(tài)儲(chǔ)氫、金屬氫化物儲(chǔ)氫等3 種儲(chǔ)氫方式各自的優(yōu)缺點(diǎn)與發(fā)展現(xiàn)狀，為未來發(fā)展提供新的思路.馬原等[33]采用CFD 方法建立兩相流模型，對(duì)微重力條件下在液氫儲(chǔ)罐內(nèi)噴射過冷流體實(shí)現(xiàn)儲(chǔ)罐降壓的過程開展數(shù)值模擬研究，對(duì)比計(jì)算了不同噴射區(qū)域、噴射流量、噴射速度等對(duì)罐內(nèi)溫度場(chǎng)分布與壓力變化的影響，認(rèn)為氣--液相區(qū)噴射降壓性能優(yōu)于單獨(dú)區(qū)域噴射，液相區(qū)噴射降壓效果最弱.2019 年王舜浩等[34]利用數(shù)值方法研究了液氫縮比儲(chǔ)罐內(nèi)部流體的蒸發(fā)性質(zhì)，通過與已有的實(shí)驗(yàn)數(shù)據(jù)對(duì)比構(gòu)建了基于VOF 兩相流模型以及(level-set)界面跟蹤方法的儲(chǔ)罐內(nèi)部流體流動(dòng)和相變傳熱傳質(zhì)模型框架，為模擬液氫儲(chǔ)罐地面停放階段的熱物理過程提供了參考.Zuo 等[35]建立了一個(gè)噴嘴可旋轉(zhuǎn)的三維液氫儲(chǔ)罐模型，利用CFD 方法研究了低溫射流對(duì)零重力條件下熱分層的影響，結(jié)果表明可旋轉(zhuǎn)的噴嘴可以很有效的抑制罐體內(nèi)部的熱分層現(xiàn)象.Guo 等[36]采用二維縮比模型儲(chǔ)罐研究了低溫射流消除熱分層現(xiàn)象中噴頭形狀、位置以及射流速度對(duì)消除效果的影響，認(rèn)為圓形射流噴嘴相較于半球形射流噴嘴消除效果更好.Zhang 等[37]研究了一種效率高,界面清晰,適用于三維模型的計(jì)算氣液兩相界面遷移特性的歐拉運(yùn)動(dòng)界面追蹤方法,該方法將‘米’狀相鄰單元Youngs 方法用于運(yùn)動(dòng)界面重構(gòu),將Youngs-VOF 和水平集通過幾何方法耦合,提高運(yùn)動(dòng)界面精度,克服了VOF 和水平集方法存在的缺陷,避免了利用高階導(dǎo)數(shù)本身的穩(wěn)定性去求解水平集對(duì)流方程和距離函數(shù)方程，為今后將熱分層現(xiàn)象的研究拓展到三維以及兩相流提供了思路.

由上可知國(guó)內(nèi)外學(xué)者已經(jīng)在該問題上進(jìn)行了很多有益嘗試并且取得了一定進(jìn)展，但是以往的研究大多假定儲(chǔ)罐壁表面均勻漏熱，因此在整個(gè)罐壁上采用均勻熱量邊界條件，并且在計(jì)算時(shí)認(rèn)為邊界漏熱和低溫射流是同時(shí)進(jìn)行的，通過對(duì)比低溫射流一段時(shí)間后儲(chǔ)罐內(nèi)部的最高溫度判斷低溫射流系統(tǒng)增壓控制性能的優(yōu)劣.儲(chǔ)罐罐體一般為金屬材料并且在表面覆蓋多層絕熱材料，盡管如此，實(shí)際使用中還是無法完全消除熱量滲入.罐體表面因?yàn)榻Y(jié)構(gòu)裝置等的差異，通過儲(chǔ)罐壁滲入液體推進(jìn)劑的熱量并不是均勻分布的.相對(duì)于罐體其他位置，通過易漏熱的區(qū)域滲入罐體的熱量對(duì)推進(jìn)劑的熱分層有著更大影響.目前對(duì)于儲(chǔ)罐局部漏熱導(dǎo)致的熱分層現(xiàn)象的低溫射流消除系統(tǒng)的設(shè)計(jì)分析研究還不夠充分，因此本文主要研究了大尺寸儲(chǔ)罐在局部區(qū)域漏熱情形下儲(chǔ)罐內(nèi)部出現(xiàn)的明顯熱分層，并對(duì)比研究了低溫射流噴嘴的形狀及其在儲(chǔ)罐內(nèi)部的相對(duì)位置等因素對(duì)利用低溫射流消除熱分層效果的影響.

1 計(jì)算模型

1.1 物理模型

熱分層控制系統(tǒng)模型由罐體結(jié)構(gòu)和入射結(jié)構(gòu)組成,如圖1 所示.低溫射流由儲(chǔ)罐內(nèi)部的射流噴嘴射出，與儲(chǔ)罐內(nèi)部流體發(fā)生混合和熱交換；另一方面儲(chǔ)罐內(nèi)部流體從罐體出流口導(dǎo)出，經(jīng)過儲(chǔ)罐外制冷系統(tǒng)(本文忽略)處理成低溫流體重新流回儲(chǔ)罐內(nèi)部，形成一個(gè)消除熱分層的閉環(huán)機(jī)制.本文主要研究罐體局部漏熱的情況，并考慮罐體出流口作為儲(chǔ)罐與外部系統(tǒng)的銜接段易發(fā)生漏熱的工程經(jīng)驗(yàn)，假定罐體表面的條狀區(qū)域及出口為漏熱帶(見圖2)，罐體其他區(qū)域視為絕熱情況.本文數(shù)值模擬采用二維軸對(duì)稱模型，儲(chǔ)罐各部分具體尺寸見圖2；采用液氫為研究工質(zhì)，填充率保持為100%.本文通過改變射流噴嘴的形狀以及射流噴嘴在儲(chǔ)罐內(nèi)部的相對(duì)位置，研究不同低溫射流條件對(duì)儲(chǔ)罐內(nèi)部流體流動(dòng)和溫度分布時(shí)空演化過程的影響.

圖1 罐體結(jié)構(gòu)及射流噴嘴樣式示意圖Fig.1 Schematic diagram of tank structure and jet nozzle pattern

圖2 儲(chǔ)罐幾何尺寸及監(jiān)測(cè)線示意圖Fig.2 Schematic diagram of the cryogenic storage tank and inspection line

1.2 數(shù)學(xué)模型

1.2.1 控制方程

假設(shè)流體是不可壓縮的并且具有恒定的熱物理特性.在微重力條件下忽略重力和浮力的影響，因而在圓柱坐標(biāo)系中整個(gè)計(jì)算域內(nèi)的質(zhì)量守恒、動(dòng)量守恒和能量守恒控制方程如下:

其中，μ為動(dòng)力學(xué)黏性系數(shù)，由工質(zhì)物理特性決定，μt為湍流黏性系數(shù)，是空間坐標(biāo)的函數(shù)，取決于流動(dòng)狀態(tài)而非物性參數(shù)，μt為湍流黏性系數(shù),是由湍動(dòng)能k及耗散率確定ε，表達(dá)式如下所示

其中?；?shù)Cμ=0.09.

(3)能量方程

其中prt≈0.85.

1.2.2 湍流模型

采用雷諾時(shí)均模型對(duì)湍流進(jìn)行模擬.為了使方程封閉，引入了新的未知量湍動(dòng)能k以及耗散率ε.Abid[38]對(duì)受限沖擊低溫射流以及受限對(duì)沖低溫射流進(jìn)行了數(shù)值模擬實(shí)驗(yàn)研究，認(rèn)為采用AB 型低雷諾數(shù)k-ε 湍流模型可以很好的對(duì)低溫射流進(jìn)行數(shù)值模擬，因此本文采用AB 低雷諾數(shù)k-ε 湍流模型，其kε 模式方程以及渦黏性如下所示

1.2.3 邊界條件以及初始條件

在入射口截面上

在對(duì)稱軸上

在罐壁漏熱帶上

在罐體其他壁面以及導(dǎo)管壁面上

初始時(shí)刻，罐內(nèi)流體溫度

出流口處邊界為自由出流條件.

工質(zhì)的物理特性根參考了Guo 等[36]的研究,選取溫度為20 K 時(shí)的各項(xiàng)參數(shù)值如表1 所示.

表1 工質(zhì)的熱物理特性參數(shù)值Table 1 Thermal physical property values of liquid hydrogen

2 數(shù)值模擬

上述控制方程的求解通過國(guó)際通用的CFD 模擬軟件Fluent 17.0 進(jìn)行計(jì)算，整個(gè)計(jì)算域通過ICEM生成網(wǎng)格，在壁面附近網(wǎng)格加密，使用SIMPLE(semi-implicit method for pressure-linked equations)算法求解，迭代過程采用絕對(duì)收斂標(biāo)準(zhǔn)進(jìn)行控制，能量項(xiàng)的殘差值設(shè)定為10?9，其他各項(xiàng)設(shè)為10?5.本文采用3 種網(wǎng)格數(shù)目分別為71 757，110 249，174 124，開展了網(wǎng)格無關(guān)性檢驗(yàn).如圖2 所示，在計(jì)算域上溫度梯度較大的位置選取了水平監(jiān)測(cè)線和豎直監(jiān)測(cè)線.圖3 顯示了網(wǎng)格無關(guān)性檢驗(yàn)的結(jié)果，可以發(fā)現(xiàn)沿豎直監(jiān)測(cè)線3 種網(wǎng)格數(shù)下的溫度分布基本一致；沿水平監(jiān)測(cè)線兩種較細(xì)網(wǎng)格的溫度分布沒有明顯差異，但在靠近低溫射流噴嘴位置處與粗網(wǎng)格下的溫度分布有明顯差異.由此，之后的研究采用網(wǎng)格數(shù)量110 426 以同時(shí)兼顧計(jì)算準(zhǔn)確性和計(jì)算量.

圖3 網(wǎng)格無關(guān)性檢驗(yàn)Fig.3 Grid independence test

3 結(jié)果討論

主要研究了圓形射流噴嘴(C)在儲(chǔ)罐內(nèi)部不同位置時(shí)利用低溫射流消除熱分層的效果，并選取了射流噴嘴所在典型位置與半球形射流噴嘴(H)進(jìn)行了對(duì)比.表2 列出了各工況進(jìn)行數(shù)值計(jì)算時(shí)的具體參數(shù)設(shè)置.

表2 各計(jì)算工況具體參數(shù)Table 2 Details of computational cases

3.1 圓形射流噴嘴

采用圓形射流噴嘴進(jìn)行研究，其中低溫射流速度沿噴嘴表面均勻分布，將無射流情況下漏熱持續(xù)30 d(2 592 000 s)后儲(chǔ)罐內(nèi)部的溫度場(chǎng)分布作為射流開始的初始狀態(tài)(t=0 s).如圖4 所示，可以發(fā)現(xiàn)此時(shí)儲(chǔ)罐內(nèi)部流體以漏熱帶為中心出現(xiàn)了明顯的環(huán)形熱分層現(xiàn)象.射流噴嘴在儲(chǔ)罐內(nèi)部的相對(duì)位置對(duì)初始熱分布沒有明顯影響，4 種工況中的熱分布基本一致.

圖4 各工況初始時(shí)刻溫度云圖分布Fig.4 Isothermals of the initial state(t=0 s)for the cryogenic jet

對(duì)于在漏熱帶施加溫度邊界條件模擬儲(chǔ)罐壁面局部漏熱的情形，漏熱帶上的熱量傳輸速率可以表征漏熱帶附近的溫度分布.熱量傳輸速率指單位時(shí)間通過漏熱帶進(jìn)入罐體的熱量，其值越大說明漏熱帶附近區(qū)域溫度越低.圖5 給出了射流時(shí)間持續(xù)約1 200 s 后各工況條狀漏熱帶熱量傳輸速率隨時(shí)間演化圖.可以看出在不同時(shí)間段，低溫射流噴嘴與漏熱帶的相對(duì)位置對(duì)于漏熱帶附近高溫區(qū)域的熱量傳輸速率存在顯著影響.按照各工況漏熱帶熱量傳輸速率的差異大體可以將射流過程分為3 個(gè)階段.第一階段(0～100 s)，對(duì)任意一種工況，由于低溫射流時(shí)間較短，通過射流進(jìn)入罐體內(nèi)部的冷流體都無法直接影響條狀漏熱帶附近區(qū)域.冷流體在罐內(nèi)作用區(qū)域非常有限，因此各工況條狀漏熱帶熱量傳輸速率曲線基本重合.圖6 給出了低溫射流進(jìn)行64.08 s 后各工況儲(chǔ)罐內(nèi)部的溫度云圖分布.可以看出低溫流體主要分布在射流噴嘴附近，距離條狀漏熱帶附近的高溫區(qū)域較遠(yuǎn).此階段儲(chǔ)罐內(nèi)抑制熱分層發(fā)展主要依靠罐內(nèi)流體由出流口經(jīng)過罐外冷卻系統(tǒng)實(shí)現(xiàn).此時(shí)出流口附近的的溫度存在明顯下降.圖7 給出了該時(shí)刻各工況的流場(chǎng)圖，可以看出施加低溫射流之后，會(huì)在射流噴嘴附近形成渦流，由于射流噴嘴位置的差異，渦流在儲(chǔ)罐內(nèi)部的位置有所不同，因此儲(chǔ)罐內(nèi)部流場(chǎng)存在差異，但總體而言在低溫射流初期，各工況抑制熱分層的效果差異不大.對(duì)比圖4 初始時(shí)刻各工況溫度云圖分布，可以發(fā)現(xiàn)條狀漏熱帶附近高溫區(qū)域的熱分層形態(tài)沒有發(fā)生明顯變化，而出口位置的熱分層明顯被消除.在第二階段(100～700 s)，各工況低溫射流持續(xù)到424.08 s 時(shí)的溫度云圖(圖8)表明工況C1 和C2 中條狀漏熱帶附近熱分層形態(tài)與第一階段相比發(fā)生了輕微的變化.由射流噴嘴進(jìn)入到罐體內(nèi)部的冷流體(溫度介于16 ～18 K 的低溫流體)作用范圍進(jìn)一步擴(kuò)大，在罐體沿著軸線形成了冷流體柱.工況C3 和C4 中條狀漏熱帶附近高溫區(qū)域的熱分層形態(tài)相較于第一階段則發(fā)生了顯著的變化.冷流體已經(jīng)在罐底聚集并沿著儲(chǔ)罐壁面內(nèi)卷.漏熱帶附近高溫區(qū)域熱分層的形態(tài)和該時(shí)刻各工況中渦流的位置有很大的關(guān)系.總體而言，熱分層的延展方向與對(duì)應(yīng)位置的流體流動(dòng)方向是一致的.由于從低溫射流開始C3 和C4 兩種工況形成的渦流相對(duì)于C1 和C2 更靠近罐體底部，長(zhǎng)時(shí)間作用之后，靠近對(duì)稱軸的熱分層隨著順時(shí)針方向流動(dòng)的流體而向左延展.圖9 給出了低溫射流進(jìn)行到424.08 s 時(shí)的流場(chǎng)圖.可以看出，各工況渦流都沿對(duì)稱軸運(yùn)動(dòng)到罐體底部.C1 中形成的渦流在罐體內(nèi)部的流動(dòng)區(qū)域范圍相對(duì)于其他工況更廣，整體性更強(qiáng)，能夠更大范圍帶動(dòng)罐內(nèi)流體的混合，結(jié)合圖5 可以發(fā)現(xiàn)這一時(shí)段，C1 中漏熱帶熱量傳輸速率高于其他工況，并且各工況整體上滿足隨著射流噴嘴伸入罐體內(nèi)部長(zhǎng)度的增加，漏熱帶熱量傳輸速率逐漸降低.這和圖10 給出的各工況儲(chǔ)罐內(nèi)部流體平均速度隨低溫射流時(shí)間演化曲線圖變化趨勢(shì)是一致的，射流初期各工況儲(chǔ)罐內(nèi)部流體平均速度都快速上升，但工況間沒有明顯差異.低溫射流持續(xù)100 s 之后，罐體內(nèi)部流體平均速度隨著射流噴嘴伸入罐體內(nèi)部長(zhǎng)度的增加而逐漸減小.圖11 給出了儲(chǔ)罐內(nèi)部流體平均溫度隨低溫射流時(shí)間的演化曲線，可以發(fā)現(xiàn)各工況儲(chǔ)罐內(nèi)部流體平均溫度的差異也主要是在低溫射流持續(xù)100 s 之后形成的，整體上各工況平均溫度都經(jīng)歷了先快速降低然后逐漸回升的過程，同樣由于儲(chǔ)罐內(nèi)部流場(chǎng)整體性的差異，在同一低溫射流時(shí)刻，隨著射流噴嘴伸入罐體內(nèi)部長(zhǎng)度的增大，罐體內(nèi)部流體平均溫度逐漸減小.

圖5 圓形射流噴嘴漏熱帶熱量傳輸速率隨時(shí)間演化圖Fig.5 Evolution diagram of the heat transfer rate with cryogenic jet time in the circular jet nozzle

圖6 各工況低溫射流64.08 s 時(shí)溫度云圖分布Fig.6 Isothermals of cases at t=64.08 s for the cryogenic jet

圖7 各工況低溫射流64.08 s 時(shí)流場(chǎng)圖Fig.7 Flow field diagram at t=64.08 s for the cryogenic jet

圖8 各工況低溫射流424.08 s 時(shí)溫度云圖分布Fig.8 Isothermals of cases at t=424.08 s for the cryogenic jet

圖9 各工況低溫射流424.08 s 時(shí)流場(chǎng)圖Fig.9 Flow field diagram at t=424.08 s for the cryogenic jet

在第三階段(700 ～1 200 s),由低溫射流持續(xù)1 024.08 s 時(shí)的溫度云圖(圖12)可以發(fā)現(xiàn)，該時(shí)刻各工況熱分層的形態(tài)都發(fā)生了顯著變化，由最初低溫射流初始時(shí)刻的環(huán)形分層演變成不規(guī)則的帶狀分層.對(duì)比第二階段的溫度云圖分布(圖8)，可以發(fā)現(xiàn)由射流噴嘴進(jìn)入罐體內(nèi)部的帶狀冷流體分布范圍并沒有繼續(xù)在罐體底部聚集或沿罐壁攀升，反而相對(duì)于上一階段有所減少，一方面是因?yàn)殡S著罐體內(nèi)部流體流速的增加(圖10)，低溫流體很快被輸運(yùn)到罐體其他位置發(fā)生熱交換，無法沿對(duì)稱軸附近積聚延伸發(fā)展，另一方面隨著漏熱時(shí)間的延長(zhǎng)以及漏熱帶熱量傳輸速率的快速升高(圖5)，罐內(nèi)流體溫度總體逐漸上升(圖11)，射流噴嘴附近的流體溫度也有所上升，通過射流噴嘴進(jìn)入罐體的冷流體在射流噴嘴附近就會(huì)發(fā)生熱交換.事實(shí)上由罐體內(nèi)部流體平均溫度的變化曲線(圖11)表明，在低溫射流持續(xù)400 s 左右時(shí)，各工況罐體內(nèi)部平均溫度都達(dá)到極小值，隨后溫度開始上升，從能量角度來看，表明在400 s 左右從儲(chǔ)罐內(nèi)部輸出的熱量和從外界環(huán)境輸入罐體內(nèi)部的熱量達(dá)到了平衡，從儲(chǔ)罐內(nèi)部輸出的熱量包括由射流噴嘴進(jìn)入罐體內(nèi)部的冷流體攜帶的熱量(負(fù)值)和通過儲(chǔ)罐出口出流的流體攜帶的熱量，而從外界環(huán)境輸入罐體內(nèi)部的熱量主要包括通過條狀漏熱帶和出口現(xiàn)階段漏熱帶進(jìn)入罐體內(nèi)部的熱量.在低溫射流初期，輸入罐體內(nèi)部的熱量值小于輸出罐體內(nèi)部的熱量，罐體內(nèi)部流體平均溫度逐漸降低，而低溫射流后期則反之.圖13 給出了各工況低溫射流1024.08 s 時(shí)流場(chǎng)圖，可以發(fā)現(xiàn)各工況形成的渦流流動(dòng)范圍都已經(jīng)擴(kuò)大到足夠覆蓋條狀漏熱帶附近的高溫區(qū)域，結(jié)合圖5 可得，在這一階段各工況條狀漏熱帶熱流傳輸速率先后開始迅速提高，這是因?yàn)楣摅w內(nèi)部的渦流演化使得渦流中心區(qū)向條狀漏熱帶附近移動(dòng)，而渦流的中心區(qū)類似于一個(gè)滯留區(qū)，從產(chǎn)生開始就攜帶著低溫流體，渦流中心區(qū)與鄰近區(qū)域的熱量交換進(jìn)行的很慢，當(dāng)其移動(dòng)靠近條狀漏熱帶邊界后，條狀漏熱帶因?yàn)楦浇鼌^(qū)域存在冷流體團(tuán)，熱量傳輸速率顯著加快.射流噴嘴越靠近罐體底部，由低溫射流引起的渦流中心區(qū)越容易到達(dá)條狀漏熱帶附近，因此各工況中條狀漏熱帶熱量傳輸速率曲線快速增長(zhǎng)的先后順序是C4,C3,C2,C1.在渦流中心區(qū)靠近條狀漏熱帶一段時(shí)間后，隨著低溫射流時(shí)間的進(jìn)行，渦流中心區(qū)溫度也會(huì)因?yàn)闊峤粨Q而上升，此后通過條狀漏熱帶進(jìn)入罐體的熱量傳輸速率就會(huì)減慢，而從罐體內(nèi)部輸出的熱量在增加，因此罐體內(nèi)部流體平均溫度就會(huì)開始下降(圖11 中C3 和C4 曲線).

圖10 圓形低溫射流噴嘴儲(chǔ)罐平均速度隨低溫射流時(shí)間演化圖Fig.10 Evolution diagram of the average velocity of the tank with the cryogenic jet time in the circular jet nozzle

圖11 圓形低溫射流噴嘴儲(chǔ)罐平均溫度隨低溫射流時(shí)間演化圖Fig.11 Evolution diagram of the average temperature of the tank with the cryogenic jet time in the circular jet nozzle

圖12 各工況低溫射流1 024.08 s 時(shí)溫度云圖分布Fig.12 Isothermals of cases at t=1 024.08 s for the cryogenic jet

圖13 各工況低溫射流1 024.08 s 時(shí)流場(chǎng)圖Fig.13 Flow field diagram at t=1 024.08 s for the cryogenic jet

通過以上分析可知4 種工況對(duì)條狀漏熱帶附近高溫區(qū)域的作用機(jī)理是一致的，但由于射流噴嘴在罐內(nèi)位置的不同，發(fā)揮作用的時(shí)間有所不同.整體而言，低溫射流抑制熱分層的機(jī)理主要有兩種，一種是在射流早期，通過促進(jìn)儲(chǔ)罐內(nèi)部流體運(yùn)動(dòng)，避免熱量在儲(chǔ)罐內(nèi)部局部區(qū)域大量積聚形成高溫區(qū)域，包括上文中第一和第二階段；另一種則是低溫射流持續(xù)發(fā)展一段時(shí)間，因?yàn)閮?chǔ)罐內(nèi)部流場(chǎng)演化，射流引起的攜帶著低溫流體的渦流中心移動(dòng)到條狀漏熱帶附近與高溫區(qū)域進(jìn)行熱量交換.結(jié)合各工況各時(shí)間段的流場(chǎng)圖，可以發(fā)現(xiàn)在射流初期，渦流中心區(qū)位于射流噴嘴附近并沿著中軸線向罐底運(yùn)動(dòng)，此時(shí)低溫射流對(duì)條狀漏熱帶附近高溫區(qū)域影響較小，然后罐體內(nèi)部流場(chǎng)逐步發(fā)展，渦流中心區(qū)沿著罐體壁面爬升向條狀漏熱帶移動(dòng)，因此條狀漏熱帶熱量傳輸速率有了顯著提高，在渦流中心區(qū)吸收從漏熱帶進(jìn)入的熱量逐漸升溫之后，漏熱帶的熱量傳輸速率會(huì)逐漸穩(wěn)定.通過對(duì)圖5 中各條曲線對(duì)時(shí)間積分得到各工況漏熱帶傳輸熱量隨時(shí)間演化圖(圖14)，可以看出工況C4 在整個(gè)低溫射流時(shí)間內(nèi)轉(zhuǎn)移熱量的能力更強(qiáng).結(jié)合圖15 給出的溫度區(qū)間面積占比累計(jì)圖(溫度區(qū)間面積占比是指將計(jì)算域溫度范圍每隔0.5 K 設(shè)置為一個(gè)溫度區(qū)間，統(tǒng)計(jì)位于各溫度區(qū)間的計(jì)算域面積占總計(jì)算域面積的百分比并繪制成累計(jì)圖)發(fā)現(xiàn)在低溫射流結(jié)束時(shí)，各工況介于21.5 ～25 K 之間的溫度占比基本一致，但是19 ～21.5 K 之間的溫度占比C4 明顯低于其他工況，也就是說明工況C4 消除高溫區(qū)域的效果更徹底，作用更明顯，這和通過漏熱帶傳輸熱量作為判據(jù)得出的結(jié)論是一致的.

圖14 圓形射流噴嘴各工況漏熱帶傳輸熱量隨時(shí)間演化圖Fig.14 Evolution diagram of heat flux over time for each circular nozzle case

圖15 圓形射流噴嘴各工況溫度區(qū)間面積占比累計(jì)圖Fig.15 Comparison of the cumulative curve of the temperature-area-ratio for each circular nozzle case

3.2 兩種射流噴嘴形狀對(duì)比

基于小尺寸儲(chǔ)罐(直徑70 mm)的研究[38]表明：圓形射流噴嘴因?yàn)槌隽鞣较蚋?，罐體內(nèi)部流場(chǎng)演化更劇烈，因此消除熱分層的效果優(yōu)于半球形射流噴嘴.在結(jié)果討論中發(fā)現(xiàn)對(duì)于大尺寸液氫儲(chǔ)罐，圓形射流噴嘴在罐內(nèi)位置對(duì)熱分層的消除機(jī)理是一致的，只是作用時(shí)間存在差異.相對(duì)而言，靠近儲(chǔ)罐出口的位置能包含其他工況罐體內(nèi)部流體流場(chǎng)演變過程，因此以兩種射流噴嘴在靠近出口的位置時(shí)消除熱分層的效果為例進(jìn)行對(duì)比分析,即工況C1 和H1.半球形射流噴嘴截面半徑同為5.5 cm，為了保持圓形射流噴嘴相同入射質(zhì)量流量，將射流速度設(shè)為0.05 m/s,其他設(shè)定與圓形射流噴嘴情況一致.計(jì)算并繪制兩種射流噴嘴漏熱帶熱量傳輸速率隨低溫射流時(shí)間演化曲線(圖16)通過對(duì)比可以發(fā)現(xiàn),圓形射流噴嘴消除熱分層的效果明顯高于半球形射流噴嘴.圖17 和圖18 分別給出半球形射流噴嘴在低溫射流進(jìn)行424.08 s 時(shí)的溫度云圖以及流場(chǎng)圖，從溫度云圖可以看出低溫流體都聚集在射流噴嘴附近，對(duì)條狀漏熱帶附近高溫區(qū)域影響很小.從流場(chǎng)圖可以發(fā)現(xiàn)半球形射流噴嘴形成的渦流發(fā)展緩慢，無法將低溫流體輸運(yùn)到條狀漏熱帶附近的高溫區(qū)域，罐體內(nèi)部熱分層抑制主要依靠流體流動(dòng)使得熱量無法聚積，因此消除熱分層的效果不如圓形射流噴嘴.

圖16 圓形與半球形射流噴嘴熱分層消除效果對(duì)比圖Fig.16 Contrast diagram of temperature stratification of circular and hemispherical jet nozzles

圖17 工況H1 在低溫射流424.08 s 的溫度云圖Fig.17 Isothermals of H1 at t=424.08 s for the cryogenic jet

圖18 工況H1 在低溫射流424.08 s 的流場(chǎng)圖Fig.18 Flow field diagram of H1 at t=424.08 s for the cryogenic jet

4 結(jié)論

本文利用軸對(duì)稱的具有低溫射流裝置的零蒸發(fā)儲(chǔ)罐模型，通過數(shù)值模擬研究了微重力條件下液氫儲(chǔ)存過程中利用低溫射流消除環(huán)境漏熱引起的熱分層現(xiàn)象的效果.通過分析儲(chǔ)罐內(nèi)部流場(chǎng)流動(dòng)和溫度分布的時(shí)間演化過程，研究了射流噴嘴在罐內(nèi)不同位置對(duì)利用低溫射流消除儲(chǔ)罐內(nèi)部熱分層效果的影響.研究結(jié)果表明：對(duì)于大尺寸儲(chǔ)罐，當(dāng)采用圓形射流噴嘴且低溫射流條件相同時(shí)，射流噴嘴的位置對(duì)罐體內(nèi)部熱分層消除效果影響不是很明顯，在本文的入射條件下，當(dāng)?shù)蜏厣淞髦脫Q率達(dá)到2%，即低溫射流時(shí)間持續(xù)700 s 時(shí)，罐體內(nèi)部熱分層的消除效果最顯著；同時(shí)，當(dāng)射流噴嘴位于儲(chǔ)罐內(nèi)部同一相對(duì)位置且入射流量相同時(shí)，圓形射流噴嘴因出流方向更集中，罐內(nèi)流場(chǎng)演變更快，熱分層消除效果比半球形射流噴嘴更好.