摘 要：為獲得對液氫地面無排氣加注性能及其影響因素的系統(tǒng)認(rèn)識(shí)，建立二維軸對稱液氫無排氣加注數(shù)值模型，借助Fluent軟件研究不同進(jìn)液溫度、進(jìn)液流量、初始壁溫和罐內(nèi)初始?jí)毫σ簹錈o排氣加注特性的影響規(guī)律。采用Lee模型和用戶自定義程序（UDF）描述液氫在加注中的蒸發(fā)、冷凝等相變問題，并考慮氣液界面張力。結(jié)果表明：在地面加注條件下，液氫儲(chǔ)罐內(nèi)膽中氣相空間的壓力在加注初期快速升高，然后趨于平穩(wěn)，到加注后期再次明顯增大。當(dāng)所加入液體溫度對應(yīng)的飽和蒸氣壓大于容器內(nèi)部壓力時(shí)，初期加注將導(dǎo)致閃蒸現(xiàn)象，初始階段壓升顯著；進(jìn)液流量不僅影響加注總時(shí)間，還影響加注初期的汽化增壓、體系溫降程度及加注末期的氣氫壓縮效應(yīng)；初始壁溫影響加注初期的壁面沸騰增壓效應(yīng)；初始?jí)毫幼⒊跗诘膲毫ψ兓绊戄^大，當(dāng)初始?jí)毫Ω哂谝簹滹柡蛪毫r(shí)，罐內(nèi)不發(fā)生閃蒸。

關(guān)鍵詞：液氫；常重力；無排氣加注；加注特性；數(shù)值模擬；氣相空間壓力

中圖分類號(hào)：TQ116.2+7 " 文獻(xiàn)標(biāo)志碼：A

0 引 言

隨著不可再生能源的持續(xù)過量消耗，全球能源和環(huán)境的可持續(xù)發(fā)展面臨越來越嚴(yán)峻的考驗(yàn)，開發(fā)利用綠色高效無污染的替代能源成為世界各國的廣泛共識(shí)和亟需解決的問題［1］。氫作為一種清潔、高效和可持續(xù)的二次能源，具有熱值高、來源廣泛、儲(chǔ)量豐富等優(yōu)點(diǎn)，且能作為中間介質(zhì)實(shí)現(xiàn)電、氣、熱和工業(yè)原料的相互轉(zhuǎn)化與耦合，被認(rèn)為是21世紀(jì)的終極能源［2-4］。高效儲(chǔ)運(yùn)技術(shù)對于氫能產(chǎn)業(yè)鏈走向成熟并保持穩(wěn)健發(fā)展態(tài)勢具有突出重要性，不可或缺，目前還有一些瓶頸問題需要解決。與高壓氣態(tài)氫氣儲(chǔ)存技術(shù)相比，低溫液氫儲(chǔ)存的壓力大大降低，而體積和質(zhì)量儲(chǔ)氫密度顯著增大，且液氫通過專用管道進(jìn)行加注或轉(zhuǎn)注時(shí)質(zhì)量流量大、用時(shí)短，適合遠(yuǎn)距離大宗運(yùn)輸。在航空航天領(lǐng)域，這些優(yōu)點(diǎn)使得低溫液態(tài)氫氣儲(chǔ)運(yùn)技術(shù)發(fā)揮出巨大作用。隨著2021年11月GB/T 40045—2021《氫能汽車用燃料 液氫》、GB/T 40060—2021《液氫貯存和運(yùn)輸技術(shù)要求》和GB/T 40061—2021《液氫生產(chǎn)系統(tǒng)技術(shù)規(guī)范》3項(xiàng)國家標(biāo)準(zhǔn)正式頒布實(shí)施［5-7］，面向民用領(lǐng)域的低溫液氫儲(chǔ)運(yùn)技術(shù)發(fā)展開始逐步提速。在未來的民用領(lǐng)域，向固定式或移動(dòng)式液氫罐中加注液氫將成為一個(gè)常規(guī)操作。為實(shí)現(xiàn)安全、高效、經(jīng)濟(jì)地加注液氫，需要對影響液氫加注特性的主要工藝條件進(jìn)行系統(tǒng)研究，掌握液氫加注工藝的優(yōu)化方法。

國內(nèi)外學(xué)者通過理論分析［8-9］、實(shí)驗(yàn)測試［10-12］和數(shù)值模擬［13-15］等方法，對液氮、液氧、液化天然氣、液氦等不同種類低溫液體的加注過程、加注工藝和加注特性進(jìn)行了較深入研究，分析了不同設(shè)備結(jié)構(gòu)、加注工藝和加注環(huán)境條件下的加注過程。Moran等［11］進(jìn)行了液氮和液氫兩種介質(zhì)的無排氣加注實(shí)驗(yàn)研究，發(fā)現(xiàn)頂部加注和底部加注所導(dǎo)致的罐內(nèi)壓力變化特征明顯不同； 王彩莉等［12］以液氮為實(shí)驗(yàn)介質(zhì)，開展了地面加注實(shí)驗(yàn)，對比分析了無排氣加注和排氣加注過程中的性能參數(shù)變化規(guī)律； 馬原等［16］在對微重力條件下液氫的無排氣加注過程進(jìn)行模擬研究時(shí)，通過用戶自定義函數(shù)（user defined function，UDF）考慮了閃蒸、蒸發(fā)和冷凝等相變過程；雙婧婧等［17］為研究液氫儲(chǔ)罐在微重力條件下減壓過程中壓力波動(dòng)和熱效率的動(dòng)態(tài)特性，基于集總蒸汽假設(shè)，建立液氫儲(chǔ)罐的瞬態(tài)CFD模型，對液氫罐內(nèi)過冷液體的混合過程進(jìn)行了研究，通過編寫UDF程序，求解了液相和氣相之間的相變與傳熱。

液氫的沸點(diǎn)在常壓下約為-253 ℃，比環(huán)境溫度低得多，從環(huán)境漏入內(nèi)膽中的熱量會(huì)造成液氫汽化，引起內(nèi)膽中氣相空間壓力升高。在液氫加注和儲(chǔ)存過程中，內(nèi)膽中累計(jì)過量的液氫汽化和壓力升高，會(huì)影響儲(chǔ)氫系統(tǒng)的壓力穩(wěn)定性和溫度穩(wěn)定性。液氫加注過程中，液態(tài)氫和氣態(tài)氫在密度、比熱容、黏度等理化性質(zhì)上對壓力和溫度等熱力學(xué)條件非常敏感［18］。這些特性導(dǎo)致迄今為止，液氫加注時(shí)內(nèi)膽中復(fù)雜的兩相流及其演化特征等關(guān)鍵問題尚不十分清楚［19］。

本文以一臺(tái)內(nèi)膽容積為34 L的小型立式液氫儲(chǔ)罐為研究對象，建立零持液率液氫罐頂部無排氣加注的二維幾何與網(wǎng)格模型，基于CFD仿真技術(shù)，對進(jìn)液溫度、進(jìn)液流量、初始壁溫和初始?jí)毫?種加注工藝參數(shù)進(jìn)行單因素模擬研究，旨在探究加注工藝對液氫無排氣加注性能的影響，以期為液氫無排氣加注工藝優(yōu)化提供指導(dǎo)和參考。

1 仿真模型及驗(yàn)證

1.1 幾何模型

圖1是液氫儲(chǔ)罐內(nèi)膽幾何模型，材質(zhì)為S30408不銹鋼，相關(guān)尺寸參考文獻(xiàn)［11］設(shè)定為內(nèi)徑318 mm，壁厚3 mm，總高度為514 mm。加注口直徑為10 mm，位于上封頭正中心處。為避免液氫離開加注口后直沖向底部，在圖示位置增加直徑為10 mm的圓形擋板，使進(jìn)入的液氫能一定程度上沿徑向噴出。本文僅對內(nèi)膽中的液氫加注過程進(jìn)行CFD仿真研究，故幾何模型只包括內(nèi)膽和頂部的加注結(jié)構(gòu)，不含加注管路、支撐結(jié)構(gòu)和安全附件等。由于液氫儲(chǔ)罐和加注口結(jié)構(gòu)的對稱性，為降低對計(jì)算資源的需求，構(gòu)建二維軸對稱幾何模型。

1.2 CFD數(shù)值模型

將仿真計(jì)算區(qū)域分為流體域和固體域，使用 Gambit 2.4.6對液氫儲(chǔ)罐幾何模型進(jìn)行網(wǎng)格劃分，均劃分為結(jié)構(gòu)化網(wǎng)格，如圖2所示。

在進(jìn)行CFD仿真計(jì)算時(shí)，假定內(nèi)膽為剛性容器、總?cè)莘e保持不變。為研究流體區(qū)域內(nèi)氣液兩相流動(dòng)和相變行為，計(jì)算時(shí)選擇VOF多相流模型來追蹤氣液兩相界面的情況，選擇k-ω模型描述罐內(nèi)的湍流行為。氣態(tài)氫密度采用理想氣體模型進(jìn)行描述。液氫的密度一般是壓力和溫度的函數(shù)，但在本文仿真計(jì)算的壓力范圍內(nèi)，液氫密度隨壓力所產(chǎn)生的變化很小，而隨溫度的變化相對較顯著，故在本文中，僅考慮溫度對液氫密度的影響，認(rèn)為液氫密度是溫度的單值函數(shù)，采用Boussinesq假設(shè)實(shí)時(shí)更新液氫密度［20-21］。其他物性參數(shù)（熱導(dǎo)率、比熱容和黏度等）均采用溫度的分段線性函數(shù)進(jìn)行計(jì)算。

在液氫加注過程中，影響罐內(nèi)壓力和溫度變化的一個(gè)重要因素是相變，包括閃蒸、蒸發(fā)和冷凝。相變過程的建模和相關(guān)計(jì)算是液氫加注數(shù)值仿真研究中需要考慮的關(guān)鍵問題之一。本文采用Lee模型描述液氫在加注中的蒸發(fā)和冷凝相變問題，閃蒸現(xiàn)象以用戶自定義函數(shù)（UDF）形式編譯嵌入到Fluent中進(jìn)行仿真計(jì)算。

在液氫加注初期，如果儲(chǔ)罐內(nèi)初始?jí)毫p0]顯著小于進(jìn)口液氫溫度對應(yīng)的飽和壓力[psat]，最先進(jìn)入的液氫便會(huì)發(fā)生閃蒸。液氫閃蒸的劇烈程度用閃蒸系數(shù)[β]來描述，如式（1）所示。當(dāng)[p0lt;psat]時(shí)，進(jìn)入罐內(nèi)的液體發(fā)生閃蒸來平衡罐內(nèi)壓差，此時(shí)，氣態(tài)氫和液氫的質(zhì)量增加速率分別如式（2）和式（3）所示［22］：

[β=0，" p0≥psat hin-hl，sathfg，" p0lt;psat] （1）

[mg=βmin] （2）

[ml=1-βmin] （3）

式中：[hin]——加注液體焓，J/（kg·K）；[hl，sat]——罐內(nèi)壓力下對應(yīng)的飽和液體焓，J/（kg·K）；[hfg]——罐內(nèi)壓力下對應(yīng)的汽化潛熱，J/（kg·K）；[mg]——?dú)怏w質(zhì)量變化率，kg/s；[min]——液氫加注質(zhì)量流量，kg/s；[ml]——液體質(zhì)量變化率，kg/s。

對于蒸發(fā)和冷凝問題，假設(shè)相變發(fā)生在恒壓和準(zhǔn)熱平衡狀態(tài)，傳質(zhì)主要取決于飽和溫度。根據(jù)Lee模型［23］通過比較單元內(nèi)溫度[T]和對應(yīng)壓力下的飽和溫度[Tsat]，判定是發(fā)生蒸發(fā)還是發(fā)生冷凝。如果[Tgt;Tsat]，液體蒸發(fā)，否則蒸汽冷凝。蒸發(fā)和冷凝過程的傳質(zhì)過程可通過式（4）進(jìn)行求解：

[Sm=rlαlρlT-TsatTsat ， T≥TsatrgαgρgTsat-TTsat， Tlt;Tsat] （4）

式中：[Sm]——?dú)庖合嘧儺a(chǎn)生的質(zhì)量源項(xiàng)，kg/（m3·s）；[r]——經(jīng)驗(yàn)參數(shù)（稱為傳質(zhì)強(qiáng)度因子），1/s；[α]——流體體積分?jǐn)?shù)；[ρ]——流體密度，kg/m3；下標(biāo)l和g——液相和氣相。

表面張力采用連續(xù)介質(zhì)表面力（continuum surface force，CSF）模型［24］，可用穿過表面的壓力躍變來表示；體積力[Fvol]作為動(dòng)量方程中的源項(xiàng)，用散度定理表示為：

[Fvol=σlgαlρlκg▽αg+αgρgκl▽αl0.5ρl+ρg] （5）

式中：[σlg]——?dú)庖褐g的界面張力，N；[κ]——液體和蒸氣的曲率，可表示為：

[κl=▽αl▽αl， κg=▽αg▽αg] （6）

本模型每個(gè)網(wǎng)格單元內(nèi)只有兩相，則[κl=-κg]，[▽αl=-▽αg]，故式（5）可簡化為：

[Fvol=σlgρκl▽αl0.5ρl+ρg] （7）

入口邊界選擇質(zhì)量流量，且每組仿真工況下的入口進(jìn)液流量為定值。采用PISO（the Pressure-Implicit with Splitting of Operators）算法來求解半隱式壓力修正方程。梯度、壓力和體積分?jǐn)?shù)分別采用PRESTO（Green-Guass Cell Based、Pressure Staggering Option）和Geo-Reconstruct schemes進(jìn)行離散。其余項(xiàng)均采用二階迎風(fēng)格式進(jìn)行離散，松弛因子和其他模型參數(shù)均采用默認(rèn)值。根據(jù)實(shí)際收斂性能，時(shí)間步長取0.001 s，能量方程的收斂準(zhǔn)則殘差控制在10-6，其他方程的控制在10-4。

1.3 CFD模型驗(yàn)證

由于液氫的低溫危害性和實(shí)驗(yàn)條件的苛刻，本研究以文獻(xiàn)中的液氫無排氣加注實(shí)驗(yàn)條件（見表1）和結(jié)果為參考，驗(yàn)證CFD模型的正確性。

圖3給出了罐內(nèi)壓力的仿真計(jì)算結(jié)果與文獻(xiàn)［11］中實(shí)驗(yàn)測試結(jié)果的對比?？梢娀谇拔乃隼碚撃Ｐ汀⒒炯僭O(shè)、網(wǎng)格模型及初邊值條件設(shè)定，所得仿真計(jì)算結(jié)果與實(shí)測結(jié)果趨勢一致，總體吻合度較好。在圖3中，加注中期仿真計(jì)算壓力小于實(shí)測值、加注后期仿真計(jì)算壓力大于實(shí)測值，主要原因如下：1）仿真計(jì)算中假定進(jìn)口液氫的進(jìn)液流量為定值，但實(shí)際加注過程中該值有變化；2）仿真計(jì)算中假定內(nèi)容器外表面絕熱，而在實(shí)際加注過程中有外界漏熱；3）仿真計(jì)算中將氫氣視為理想氣體，實(shí)際加注過程中，氫氣并非理想氣體；4）在模擬加注時(shí)液氫從頂部進(jìn)入后沿徑向噴淋到罐內(nèi)，與文獻(xiàn)［11］的實(shí)驗(yàn)中噴灑加注存在差異。盡管仿真計(jì)算結(jié)果與實(shí)驗(yàn)測試結(jié)果存在這些差別，但總體而言可認(rèn)為本文的仿真模型能正確描述液氫加注過程中的熱質(zhì)傳遞過程，可據(jù)此開展系統(tǒng)的仿真研究。

2 仿真內(nèi)容

根據(jù)液氫實(shí)驗(yàn)研究相關(guān)文獻(xiàn)可知，用不同容積的液氫儲(chǔ)罐開展液氫加注相關(guān)實(shí)驗(yàn)研究時(shí)，進(jìn)液溫度多設(shè)定在17～22 K范圍內(nèi)；對于進(jìn)液流量的設(shè)定，加注時(shí)間足夠短很重要，因此液氫儲(chǔ)罐容積越大，允許的進(jìn)液流量上限越大。初始壁溫與進(jìn)液溫度之差一般被限制在40～100 K之間。進(jìn)液溫度對應(yīng)的飽和壓力與初始?jí)毫χ羁刂圃?0～100 kPa 范圍內(nèi)。本文確定了16組加注工藝組合，如表2所示。

基于前述經(jīng)過驗(yàn)證的仿真計(jì)算模型開展4因素4水平仿真研究。本文對于仿真結(jié)果的分析從兩個(gè)方面展開：一是頂部無排氣充裝的一般特征，二是不同加注工藝對加注期間內(nèi)膽中壓力變化的具體影響。

3 仿真結(jié)果分析

3.1 液氫頂部無排氣加注的一般特征

圖4展示了按照表1所列參數(shù)進(jìn)行液氫充裝時(shí)，罐內(nèi)壓力隨時(shí)間變化趨勢。根據(jù)圖4的曲線變化趨勢，可將低溫液氫的無排氣加注過程分為3個(gè)階段：初始快速升壓階段（持續(xù)時(shí)間為10 s）、相對穩(wěn)定加注階段（持續(xù)時(shí)間為100 s）和末期快速升壓階段（持續(xù)時(shí)間為55 s）。在加注初期，一方面，如果罐內(nèi)壓力低于進(jìn)液溫度對應(yīng)的飽和壓力會(huì)發(fā)生閃蒸；另一方面，新注入的液氫與溫度相對較高的“熱壁面”接觸而沸騰汽化，這兩方面原因共同導(dǎo)致罐內(nèi)壓力快速上升。隨著加注過程的持續(xù)，罐內(nèi)壓力不斷升高，閃蒸效應(yīng)迅速減弱，當(dāng)罐內(nèi)壓力達(dá)到液氫溫度對應(yīng)的飽和壓力時(shí)，壓力響應(yīng)過渡到第2階段。在相對穩(wěn)定加注階段，由于過熱內(nèi)膽壁面的加熱作用，仍可能發(fā)生局部的微量蒸發(fā)，傾向于引起壓力增大。隨著加注過程持續(xù)，液氫的淹沒高度不斷增大，內(nèi)膽壁及其內(nèi)部空間持續(xù)趨冷，使得內(nèi)膽中的部分氣態(tài)氫有機(jī)會(huì)重新冷凝成液滴，引起壓力降低。如此看來，在平穩(wěn)加注階段的不同時(shí)刻，內(nèi)膽中氣相空間的壓力是升高還是降低，取決于液氫蒸發(fā)和氣態(tài)氫冷凝兩種效應(yīng)孰強(qiáng)孰弱。如果這兩種效應(yīng)引起的壓力變化接近互相抵消時(shí)，就會(huì)以某種恒壓加注的特性表現(xiàn)。隨著液氫加注量的增多，當(dāng)罐內(nèi)溫度與進(jìn)液溫度相同時(shí)，蒸發(fā)和冷凝效應(yīng)均在前一階段的基礎(chǔ)上顯著減弱，氣相空間壓力的變化轉(zhuǎn)變?yōu)橹饕芤好嫔弋a(chǎn)生的機(jī)械壓縮效應(yīng)的影響，結(jié)果是壓力再次顯著增大。

圖5是加注過程中不同時(shí)刻液氫相分布云圖。由圖5可見，隨著加注的進(jìn)行，液相空間逐漸增大，氣液界面持續(xù)升高且漸趨穩(wěn)定。在不同時(shí)刻，罐內(nèi)的瞬時(shí)相分布表現(xiàn)出不同的特點(diǎn)。加注到20 s時(shí)，液氫體積分?jǐn)?shù)為12.88%，氣相空間范圍大且溫度高于液氫的溫度，液氫在下落過程中受周圍氣態(tài)氫的對流加熱和內(nèi)膽壁面輻射加熱的聯(lián)合作用，部分汽化，導(dǎo)致下落中的“液氫柱”變得細(xì)?。灰簹渎淙氲撞康囊撼睾?，激散了正下方的部分液氫，在氣液界面產(chǎn)生波紋，而液氫下落的動(dòng)能最終以熱能形式耗散掉。當(dāng)加注到60 s時(shí)，體系整體溫度進(jìn)一步降低，液氫在下落過程中的蒸發(fā)減弱，下落中的“液氫柱”明顯較20 s時(shí)的粗，且由于液池體積明顯增大（液氫的體積分?jǐn)?shù)為38.76%），下落中的液氫對氣液界面的擾動(dòng)顯著減弱，氣液界面波動(dòng)變小。此外，液池中的藍(lán)色點(diǎn)狀微區(qū)（可從電子版查看）顯示該處有氣態(tài)氫產(chǎn)生，原因可能是液氫下落過程中表面蒸發(fā)產(chǎn)生的氫氣未來得及逃逸而被帶入液池。當(dāng)加注時(shí)間為100 s時(shí)，氣液界面進(jìn)一步趨于穩(wěn)定，液池內(nèi)的二次相變現(xiàn)象也進(jìn)一步減弱。

圖6是不同加注時(shí)刻內(nèi)膽中的溫度分布云圖。由圖6可見，罐內(nèi)溫度形成了一定的分層，但隨著加注過程的持續(xù)，分層的程度逐漸減弱。溫度分層的減弱體現(xiàn)在兩個(gè)方面：一是最大溫差逐漸減小，由20 s時(shí)的40.3 K逐漸減小到160 s時(shí)的29.8 K；二是相對恒溫的區(qū)域（液相區(qū)）逐漸增大，壓縮了溫度分層區(qū)域的空間。從溫度分層的特點(diǎn)來看，不同時(shí)刻的溫度云圖大致均可細(xì)分成3個(gè)子區(qū)，自下而上分別為液相區(qū)、中間氣相區(qū)和頂部氣相區(qū)。其中中間氣相區(qū)的溫差大致在10 K以內(nèi)，頂部氣相區(qū)溫差較大且分層明顯，但區(qū)域占比較小。60 s之后，液面升高對氣相空間的擠壓作用顯著增強(qiáng)，對流作用有所加劇，這從60 s之后的溫度云圖中氣相空間有明顯的非水平等溫線可看出。加注結(jié)束時(shí)，罐內(nèi)被液氫占據(jù)的空間溫度和液氫一致，只有頂部氣相空間的溫度高于液氫溫度。從安全角度講，較大的溫度分層會(huì)對液氫加注產(chǎn)生不利影響，為減弱溫度分層程度，需要在液氫加注前進(jìn)行充分預(yù)冷，還應(yīng)盡量避免采用長徑比大的立式儲(chǔ)罐。

3.2 進(jìn)液溫度對加注性能的影響

在液氫加注過程中，進(jìn)液溫度是一個(gè)重要的工藝參數(shù)。圖7給出了進(jìn)液溫度分別為18、19、20和21 K時(shí)罐內(nèi)壓力隨充裝過程的變化情況，其他條件設(shè)置如表2（[T1、T2、T3、T4]）所示。在圖7中，由于不同加注條件對應(yīng)的加注時(shí)間有變化，故不同加注曲線的結(jié)束點(diǎn)有所不同。4種進(jìn)液溫度對應(yīng)的飽和壓力（絕壓）分別為46.6、66.0、90.7和121.5 kPa，而罐內(nèi)初始?jí)毫?6 kPa，顯著低于4種進(jìn)液溫度下對應(yīng)的液氫飽和壓力，意味著這4種情況下加注初始階段均會(huì)發(fā)生閃蒸。

由圖7可知，罐內(nèi)初始?jí)毫σ欢〞r(shí)，進(jìn)液溫度越高，閃蒸階段壓力變化越快，原因是進(jìn)液溫度越高，閃蒸系數(shù)越大，閃蒸效應(yīng)越強(qiáng)烈。由圖7還可看出：進(jìn)液溫度越高，相對穩(wěn)定加注階段的持續(xù)時(shí)間越短。當(dāng)進(jìn)液溫度從18 K升高到20 K時(shí)，穩(wěn)定加注時(shí)間縮短了57%；當(dāng)進(jìn)液溫度為21 K時(shí)，壓力變化曲線明顯和另外3條不同，即壓力坪臺(tái)區(qū)消失。產(chǎn)生這種現(xiàn)象的原因是：由于進(jìn)液溫度相對較高，液氫下落過程中引起的降溫效應(yīng)相對較弱，導(dǎo)致加注中期的蒸發(fā)-冷凝競爭格局難以達(dá)到均衡狀態(tài)，即蒸發(fā)效應(yīng)始終強(qiáng)于冷凝效應(yīng)，因此壓力曲線持續(xù)走高。在加注后期，氣相溫度進(jìn)一步降低、壓力因機(jī)械壓縮效應(yīng)占主導(dǎo)地位而再次快速升高。當(dāng)氣相溫度偏低而壓力偏高時(shí)，壓力預(yù)測結(jié)果將明顯偏大，因?yàn)闅怏w性質(zhì)已明顯偏離理想氣體模型了。

本文根據(jù)美國國家標(biāo)準(zhǔn)與技術(shù)研究院（National Institute of Standards and Technology，NIST）的物性計(jì)算軟件得到的氣態(tài)氫的壓縮因子判斷理想氣體模型的適用范圍。根據(jù)NIST的物性計(jì)算軟件，當(dāng)氣氫的壓力和溫度分別為400 kPa和28 K時(shí)，壓縮因子為0.8182；當(dāng)氣氫的壓力和溫度分別為600 kPa和32 K時(shí)，壓縮因子為0.8076。這意味著在該參數(shù)范圍內(nèi)，采用理想氣體狀態(tài)模型引起的誤差將達(dá)到20%。為保證仿真計(jì)算結(jié)果的合理性和可信度，本文在處理數(shù)據(jù)時(shí)舍棄模擬計(jì)算結(jié)果中理想氣體模型引起的誤差超過20%的部分。事實(shí)上，在加注初始和近似恒壓加注階段（氣相空間壓力在200 kPa以內(nèi)），氣相空間采用理想氣體狀態(tài)方程引起的誤差在10%以內(nèi)。

3.3 進(jìn)液流量對加注性能的影響

圖8給出進(jìn)液流量分別為0.025、0.035、0.045和0.065 kg/s時(shí)內(nèi)膽中壓力隨時(shí)間的變化情況，其他條件設(shè)置如表2（[T5、T6、T7、T8]）所示。由圖8可知，4條壓力變化曲線除在時(shí)間跨度方面有明顯差異外，變化趨勢基本一致。進(jìn)液流量從0.025 kg/s增加到0.065 kg/s時(shí)，加注總時(shí)間縮短61%。

在實(shí)際液氫加注中，雖然加注時(shí)間短是首要和基本的要求，但在其他條件給定時(shí)增大進(jìn)液流量，一方面會(huì)使加注初期的汽化更劇烈和壓升更迅猛，另一方面也會(huì)在一定程度上增強(qiáng)加注后期的氣體壓縮效應(yīng)，影響加注終了的氣相區(qū)壓力和溫度。若在未達(dá)到額定充滿率之前，罐內(nèi)壓力已達(dá)到了設(shè)定壓力，則充裝將不能繼續(xù)，從而造成充液率不達(dá)標(biāo)的情況，因此進(jìn)液流量不宜過大。

3.4 初始壁溫對加注性能的影響

圖9給出了初始壁溫分別為30、70、90和120 K時(shí)加注過程中罐內(nèi)壓力隨時(shí)間的變化曲線。其他條件設(shè)置如表2（[T9、T10、T11、T12]）所示。由圖9可知，在4種初始壁溫下，壓力變化趨勢一致，但罐內(nèi)壓力峰值隨著初始壁溫的升高而增大，初始壁溫從30 K增加到120 K時(shí)，壓力峰值顯著提高。這主要是因?yàn)槌跏急跍卦礁撸薇趦?chǔ)存的熱量越多，在冷卻時(shí)需要更多的蒸發(fā)來吸收熱量，導(dǎo)致壁面沸騰蒸發(fā)生成的氣體更多，壓力升高更顯著。可以合理地推測，如果初始壁溫過高，早期的沸騰會(huì)直接將壓力提高到上限，導(dǎo)致無排氣加注失敗。故在給定條件下，存在一個(gè)臨界初始壁溫，溫度低于該臨界值時(shí)才能實(shí)現(xiàn)有效加注。馬原等［16］在進(jìn)行液氫底部無排氣加注時(shí)也得到了同樣的結(jié)論。

3.5 初始?jí)毫幼⑿阅艿挠绊?/p>

圖10給出了罐內(nèi)初始?jí)毫Γń^壓）分別為30、40、60和80 kPa時(shí)加注過程中罐內(nèi)壓力隨時(shí)間的變化情況，其他條件設(shè)置如表2（[T13、T14、T15、T16]）所示。由圖10可知：在加注前期，不同的初始?jí)毫迌?nèi)壓力變化影響較大，當(dāng)初始?jí)毫?0 kPa時(shí)，大于進(jìn)液溫度對應(yīng)的飽和壓力（66 kPa），故在加注前期不發(fā)生閃蒸效應(yīng)，即無閃蒸加注。在加注過程中落在內(nèi)膽壁面的液氫會(huì)在熱壁面的加熱下沸騰汽化，傾向于引起壓力升高，但因液氫加注引起的體系溫度降低所導(dǎo)致的氣相壓力下降程度更大，所以在加注初期壓力不升反降。當(dāng)壓力下降到66 kPa時(shí)（液氫飽和壓力），壓力曲線過渡到第2階段。

雖然初始?jí)毫σ簹浼幼⒅泻笃趬毫蜏囟扔绊懖伙@著，且增大初始?jí)毫τ兄跍p小初始?jí)翰睿瑥亩鴾p弱閃蒸效應(yīng)，但罐內(nèi)過高的初始?jí)毫ο喈?dāng)于增大了供氫系統(tǒng)的背壓，帶來增壓加注系統(tǒng)能耗的增大。因此，初始?jí)毫σ膊灰诉^大。

4 結(jié) 論

本文以一臺(tái)內(nèi)膽水容積為34 L的小型立式液氫儲(chǔ)罐為研究對象，基于CFD數(shù)值模擬技術(shù)，對進(jìn)液溫度、進(jìn)液流量、初始壁溫和初始?jí)毫?個(gè)加注工藝條件進(jìn)行單因素模擬研究，主要結(jié)論如下：

1）向初始持液率為0（空罐）的液氫儲(chǔ)罐中加注液氫時(shí)，內(nèi)膽中氣相空間的壓力明顯依次經(jīng)歷快速升高、平穩(wěn)變化和再次攀升3個(gè)階段。加注初期壓力快速升高的原因是液氫減壓閃蒸和過熱壁面引起液氫蒸發(fā)；在相對穩(wěn)定加注階段，液態(tài)氫的蒸發(fā)以及氣態(tài)氫的再次冷凝和液面上升所致的壓縮效應(yīng)，是容器內(nèi)壓力和熱質(zhì)傳遞過程中動(dòng)態(tài)變化的耦合影響因素。在加注后期，液面的上升導(dǎo)致對氣相空間的擠壓作用成為主要影響壓力變化的因素。

2）進(jìn)液溫度、進(jìn)液流量、初始壁溫和初始?jí)毫幼⑿阅苡胁煌挠绊?。進(jìn)液溫度越高，則加注前期閃蒸越劇烈，初始階段的壓升和溫降效應(yīng)越強(qiáng)；進(jìn)液流量直接影響加注完成的時(shí)間（當(dāng)進(jìn)液流量從0.025 kg/s增加到0.065 kg/s時(shí)，加注時(shí)間從90 s減少到35 s），還影響加注初始階段的閃蒸增壓效應(yīng)及體系溫降程度和加注末期的氣氫壓縮效應(yīng)；初始壁溫主要影響加注初期的壁面沸騰增壓效應(yīng)和閃蒸末期的壓力峰值，初始壁溫過高，則早期的液氫沸騰會(huì)使壓力很快達(dá)到壓力上限，造成加注失??；初始?jí)毫χ饕绊懠幼⑶捌诘膲毫ψ兓ㄅc進(jìn)液溫度對應(yīng)的飽和蒸氣壓有關(guān)），當(dāng)初始?jí)毫Ω哂谝簹鋵?yīng)的飽和壓力時(shí)，罐內(nèi)不發(fā)生閃蒸，加注時(shí)罐內(nèi)壓力會(huì)先減小至飽和壓力，隨后再經(jīng)歷平臺(tái)區(qū)和末尾的上升區(qū)。

［參考文獻(xiàn)］

［1］ 耿銀良， 應(yīng)強(qiáng)， 李建立， 等. 液氫泄漏擴(kuò)散及安全預(yù)防研究進(jìn)展［J］. 真空與低溫， 2023， 29（2）： 153-162.

GENG Y L， YING Q， LI J L， et al. Research progress of liquid hydrogen leak diffusion and safety prevention［J］. Vacuum and cryogenics， 2023， 29（2）： 153-162.

［2］ 凌文， 劉瑋， 李育磊， 等. 中國氫能基礎(chǔ)設(shè)施產(chǎn)業(yè)發(fā)展戰(zhàn)略研究［J］. 中國工程科學(xué)， 2019， 21（3）： 76-83.

LING W， LIU W， LI Y L， et al. Development strategy of hydrogen infrastructure industry in China［J］. Strategic study of CAE， 2019， 21（3）： 76-83.

［3］ 陳宇， 張小玉， 張榮沛. 中國氫能產(chǎn)業(yè)鏈現(xiàn)狀及前景展望［J］. 新型工業(yè)化， 2021， 11（4）： 176-180， 182.

CHEN Y， ZHANG X Y， ZHANG R P. Present situation and prospect of hydrogen energy industry chain in China［J］. The journal of new industrialization， 2021， 11（4）： 176-180， 182.

［4］ SUN Z Y. 11 - Hydrogen energy ［M］. DUTTA S， MUSTANSAR HUSSAIN C. Sustainable fuel technologies Handbook. New York： Academic Press， 2021： 339-365.

［5］ GB/T 40045—2021， 氫能汽車用燃料 液氫［S］.

GB/T 40045—2021，F(xiàn)uel specification for hydrogen powered vehicles： Liquid hydrogen（LH2）［S］.

［6］ GB/T 40060—2021， 液氫貯存和運(yùn)輸技術(shù)要求［S］.

GB/T 40060—2021， Technical requirements for storage and transportation of liquid hydrogen［S］.

［7］ GB/T 40061—2021， 液氫生產(chǎn)系統(tǒng)技術(shù)規(guī)范［S］.

GB/T 40061—2021， Technical specification for liquid hydrogen production system［S］.

［8］ HEDAYAT A， CARTAGENA W， MAJUMDAR A K， et al. Modeling and analysis of chill and fill processes for the cryogenic storage and transfer engineering development unit tank［J］. Cryogenics， 2016， 74： 106-112.

［9］ KEEFER K A， HARTWIG J. Development and validation of an analytical charge-hold-vent model for cryogenic tank chilldown［J］. International journal of heat and mass transfer， 2016， 101： 175-189.

［10］ CHATO D. Ground testing of the nonvented fill method of orbital propellant transfer-results of initial test series［C］//27th Joint Propulsion Conference. Sacramento， CA， USA， 1991： 2326.

［11］ MORAN M E， NYLAND T W， DRISCOLL S L. Hydrogen no-vent fill testing in a 34 liter （1.2 cubic foot） tank［M］//FAST RW. Advances in Cryogenic Engineering. Boston， MA： Springer， 1992： 1257-1264.

［12］ WANG C L， LI Y， DENG D， et al. Performance model of the top filling configurations for no-vent fills［J］. Journal of thermophysics and heat transfer， 2011， 25（1）： 140-146.

［13］ MELIDEO D， BARALDI D. CFD analysis of fast filling strategies for hydrogen tanks and their effects on key-parameters［J］. International journal of hydrogen energy， 2015， 40（1）： 735-745.

［14］ WANG C L， LI Y， WANG R S. Performance comparison between" "no-vent" "and" "vented" fills" in" vertical" thermal-insulated" "cryogenic" cylinders［J］." Experimental" thermal and fluid science， 2011， 35（2）： 311-318.

［15］ JIANG Y B， YU Y S， WANG Z， et al. CFD simulation of heat transfer and phase change characteristics of the cryogenic liquid hydrogen tank under microgravity conditions［J］. International journal of hydrogen energy， 2023， 48（19）： 7026-7037.

［16］ MA Y， LI Y Z， ZHU K， et al. Investigation on no-vent filling process of liquid hydrogen tank under microgravity condition［J］. International journal of hydrogen energy， 2017， 42（12）： 8264-8277.

［17］ SHUANG J J， LIU Y W. Efficiency analysis of depressurization process and pressure control strategies for liquid" " hydrogen" "storage" "system" "in" " microgravity［J］. International journal of hydrogen energy， 2019， 44（30）： 15949-15961.

［18］ 劉柏文， 徐元元， 雷剛， 等. 典型低溫推進(jìn)劑的熱力學(xué)性能參數(shù)評估［J］. 火箭推進(jìn)， 2023， 49（1）： 44-53.

LIU B W， XU Y Y， LEI G， et al. Evaluation of thermodynamic performance parameters for typical cryogenic propellant［J］. Journal of rocket propulsion， 2023， 49（1）： 44-53.

［19］ 王朝， 鄒宏偉， 陳甲楠， 等. 液氮、液氫加注系統(tǒng)數(shù)值仿真與試驗(yàn)研究［J］. 太陽能學(xué)報(bào)， 2022， 43（11）： 460-465.

WANG C， ZOU H W， CHEN J N， et al. Numerical simulation and experimental study of liquid nitrogen and liquid" hydrogen" filling" system［J］." Acta" energiae" solaris sinica， 2022， 43（11）： 460-465.

［20］ 朱宇豪， 卜玉， 劉瑞敏， 等. 車載液氫瓶瓶內(nèi)自增壓過程熱力耦合特性研究［J］. 低溫工程， 2023（5）： 61-68.

ZHU Y H， BU Y， LIU R M， et al. Study of thermodynamic coupling characteristics in self-pressurization process within liquid hydrogen tank for vehicles［J］. Cryogenics， 2023（5）： 61-68.

［21］ WAN C C， ZHU S L， SHI C Y， et al. Numerical simulation on pressure evolution process of liquid hydrogen storage tank with active cryogenic cooling［J］. International journal of refrigeration， 2023， 150： 47-58.

［22］ 馬原， 王磊， 孫培杰， 等. 微重力環(huán)境低溫流體無排氣加注過程數(shù)值研究［J］. 低溫工程， 2017（4）： 7-13， 20.

MA Y， WANG L， SUN P J， et al. Numerical investigation of no-vent fill process of cryogens under microgravity condition［J］. Cryogenics， 2017（4）： 7-13， 20.

［23］ LEE W H. A pressure iteration scheme for two-phase flow" " modeling［M］. Washington， DC： Hemisphere Publishing， 1980.

［24］ BRACKBILL J U， KOTHE D B， ZEMACH C. A continuum" "method" "for" " modeling" "surface" " tension［J］. Journal of computational physics， 1992， 100（2）： 335-354.

SIMULATION STUDY OF UNVENTED FILLING CHARACTERISTICS OF LIQUID HYDROGEN

Li Jianli1，Geng Yinliang2，Li Jingfa1，Geng Jinliang3，Yu Bo1，Wu Xiaohua1

（1. School of Mechanical Engineering， Beijing Institute of Petrochemical Technology， Beijing 102617， China；

2. Technical Center of China Classification Society Quality Certification Co.， Ltd.， Beijing 100006， China；

3. Oil and Gas Regulation and Control Center of National Petroleum and Natural Gas Pipeline Network Group Co.， Ltd.， Beijing 100013， China）

Abstract：In order to obtain a systematic understandings of the unvented filling characteristics and the influencing factors of liquid hydrogen on the ground， a two-dimensional axisymmetric numerical model was established. The influences of different inlet temperatures， inlet flow rates， initial wall temperatures， and initial pressures inside the container on the unvented filling characteristics of liquid hydrogen were investigated by using Fluent software. The Lee model and a user-defined function （UDF） were used to describe the phase change issues （evaporation， and condensation） during the unvented filling of liquid hydrogen considering the interfacial tension between gas and liquid. The results show that under normal the gravity conditions， the unvented filling process of liquid hydrogen can be divided into three phases： an initial rapid pressure rise phase， a relatively stable filling phase， and a final rapid pressure rise phase. When the saturated vapor pressure corresponding to the temperature of the added liquid is greater than the pressure inside the vessel， the initial filling will lead to flash evaporation， and the pressure rise in the initial phase is significant. The inlet flow rate affects not only the total filling time， but also the vaporization pressurization at the initial stage of unvented filling， the temperature reduction degree in the system and the compression effect of gaseous hydrogen at the later stage of filling. The initial wall temperature affects the wall boiling pressurization effect during the initial filling stage. The initial pressure has a significant impact on the pressure change during the initial filling period， and the flash evaporation does not occur in the container when the initial pressure is higher than the saturation pressure of liquid hydrogen.

Keywords：liquid hydrogen； constant gravity； unvented filling； filling characteristics； numerical simulation； gas phase space pressure