水升華器升華模式理論分析與數(shù)值仿真

廖俊元楊春信楊涵

(北京航空航天大學(xué)航空科學(xué)與工程學(xué)院，北京 100191)

1 引言

隨著中國航天技術(shù)的發(fā)展及載人登月工程的逐步推進，航天員的艙外活動需求逐漸提高。艙外航天服的熱控系統(tǒng)需要保證航天員在出艙活動時處于安全、舒適的熱環(huán)境之中。水升華器是一種利用水的升華潛熱進行散熱的消耗型相變散熱裝置，相比于輻射器等熱控器件具有體積小、質(zhì)量輕、效率高、微重力下工作可靠等優(yōu)點，是一種性能優(yōu)良的空間熱沉，十分適用于應(yīng)對小散熱面、短時大功耗任務(wù)的艙外航天服的熱控系統(tǒng)。

國際對水升華器的研究開展較早，從上世紀(jì)60年代初已有關(guān)于水升華器作為熱沉應(yīng)用于航天器熱控系統(tǒng)的相關(guān)研究，并制作了如LM-107、LM-209等一系列原理樣機進行了初步的實驗探究，如圖1所示。在上世紀(jì)70年代初至上世紀(jì)末的期間，水升華器處于大量工程應(yīng)用的階段，在首次成功應(yīng)用于美國的阿波羅探月工程中后，還在蘇聯(lián)(俄羅斯)的Orlan航天服等航天工程中得到應(yīng)用。進入新世紀(jì)后，為了滿足更高的航天工程需求，水升華器出現(xiàn)了追求防污染、良好動態(tài)性能的發(fā)展趨勢，對此NASA的研究人員提出了X-38、CIS、SDC、ISDC 4種水升華器的新構(gòu)型，如圖2所示。

圖1 早期水升華器Fig.1 Early application of sublimator

圖2 新型水升華器Fig.2 New type of sublimator

中國對水升華器的相關(guān)研究起步較晚。吳志強等最早對水升華器進行了散熱性能的理論分析及一系列實驗研究；李森等建立了多孔板單孔內(nèi)水相變過程的物理數(shù)學(xué)模型，并進行了周期模式的數(shù)值模擬；王玉瑩等采用一維熱阻網(wǎng)格模型，通過仿真得到了周期模式下界面位置、水升華器溫度等參數(shù)的變化，并對不同加熱邊界的水升華器工作性能進行了實驗研究。航天工程方面，中國關(guān)于水升華器的首次成功應(yīng)用為2008年神舟七號飛船出艙活動的“飛天”艙外航天服。

由于水的蒸發(fā)、升華等相變過程發(fā)生在水升華器內(nèi)部的給水腔、多孔板微孔等尺寸微小的區(qū)域，實驗僅能測得外部溫度等宏觀數(shù)據(jù)，因此對于水升華器工作原理的研究主要通過理論分析、數(shù)值仿真來進行。但是，目前國內(nèi)對水升華器的仿真分析較少且主要為一維數(shù)值仿真，沒有考慮水升華器不同工作模式的轉(zhuǎn)換。為了探究水升華器工作模式轉(zhuǎn)換的規(guī)律，分析不同設(shè)計參數(shù)對水升華器工作特性的影響，本文在總結(jié)前人理論分析的基礎(chǔ)上，對水升華器工作模式轉(zhuǎn)換特性進行量綱分析，通過Fluent數(shù)值仿真的方法對水升華器在低熱載荷下的穩(wěn)態(tài)工作模式(升華模式)開展研究，分析不同多孔板結(jié)構(gòu)參數(shù)、熱載荷大小下水升華器工作特性的變化規(guī)律。

2 水升華器工作過程分析

2.1 水升華器工作原理

本文研究對象為平板型水升華器，其基本結(jié)構(gòu)由供水入口、加熱板、給水腔、多孔板組成，如圖3所示。

圖3 升華模式下的水升華器Fig.3 Water sublimator working in sublimation mode

系統(tǒng)工作過程為:工質(zhì)水在一定的供水壓力下被送往高真空度的給水腔，水的壓力驟降至其三相點以下而迅速蒸發(fā)。蒸發(fā)潛熱帶走大量熱量，若此時水溫降至三相點溫度以下，則會凝固形成固態(tài)冰。外表面的冰層暴露于高真空環(huán)境下，熱載荷將通過升華潛熱排向外太空。

Hamilton Standard公司根據(jù)工質(zhì)水的不同消耗形式，將水升華器工作模式劃分為4種:升華模式、蒸發(fā)模式、混合模式、周期模式。當(dāng)熱載荷較小時，水升華器工作在升華模式。由于水蒸氣流量較小，升華界面對應(yīng)的飽和蒸氣壓低于水的三相點壓力。根據(jù)圖4，此時水不會進入多孔板的毛細孔，而是首先在多孔板與給水腔的界面上形成一定厚度的冰層(圖3)。此時，供水流量等于升華流量，升華帶走的熱量等于底面加熱及供水熱量之和，水升華器內(nèi)部達到整體的熱質(zhì)動態(tài)平衡。升華模式是水升華器最期望的工作模式[12]。

圖4 水的三相圖Fig.4 Three-phase diagram of water

當(dāng)熱載荷繼續(xù)增大，水將進入多孔板形成一定厚度的冰層。冰層在升華、熔化的共同作用下逐漸耗盡后，水繼續(xù)流入并凝固成冰層，形成了周期性的循環(huán)工作過程，即周期模式。若多孔板溫度較高，則冰層不會出現(xiàn)，工質(zhì)將以液態(tài)水蒸發(fā)的形式進行散熱，即蒸發(fā)模式。水升華器實際工作中可能會出現(xiàn)以上3種模式共存的情況，這種模式被稱為混合模式。

2.2 升華模式分析

吳志強等提出升華模式下水升華器整體處于穩(wěn)態(tài)工作的狀態(tài)，可以通過理論分析計算升華熱流、加熱面溫度等參數(shù)。假設(shè)如下:

1)多孔板、加熱板及工質(zhì)水的物性為常數(shù)；

2)不考慮結(jié)構(gòu)漏熱、輻射漏熱等熱損失；

3)多孔板毛細孔簡化為眾多同徑的豎圓孔；

4)忽略給水腔內(nèi)部水對流造成的換熱；

5)加熱板、給水腔、多孔板具有相同的面積。

升華模式下水升華器內(nèi)部的熱質(zhì)平衡如圖3中展示，圖中紅色虛線箭頭代表熱量的流動，藍色實線箭頭代表質(zhì)量的流動。定義加熱板熱流密度如式(1):

式中，

q

為底面加熱熱流，W/m；

Q

為水升華器熱載荷，W；

A

為加熱板面積，m。進一步對液態(tài)水、冰層分別進行熱平衡分析，可得式(2)、(3):

式中，

q

、

q

分別為水層、冰層的導(dǎo)熱熱流，W/m；

m

為供水流量，kg/s，在升華模式下等于總體升華流量；

t

、

t

分別為供水溫度、水-冰相變溫度，K；

c

為水的比熱，J/(kg·K)；Δ

h

、Δ

h

分別為水-冰相變潛熱及冰-水蒸氣升華潛熱，J/kg。同時，還可以分別寫出水層與冰層的導(dǎo)熱關(guān)系如式(4)、(5):

式中，

λ

、

λ

分別為冰和水的導(dǎo)熱系數(shù)，W/(m·K)；

t

、

t

為升華溫度、底面溫度，K。需要說明的是，式(4)僅考慮了水層的穩(wěn)態(tài)導(dǎo)熱，而忽略了入口較高溫度供水的影響。王玉瑩在相關(guān)仿真中做了類似的處理，其仿真與實驗結(jié)果較為一致，可以認為該處理方式的誤差在可接受范圍內(nèi)。

此外根據(jù)式(6)所示，升華模式下水層、冰層的厚度之和與給水腔厚度相等:

式中，

d

、

d

、

d

分別為冰層厚度、水層厚度、給水腔厚度，m。根據(jù)式(3)，可推導(dǎo)出水升華器消耗工質(zhì)水的質(zhì)量流量如式(7):

式中，

ρ

為水的密度，kg/m；，

S

、

u

分別為供水入口的面積與供水水流速度。由于升華質(zhì)量流量很小，水蒸氣分子平均自由程遠大于流動特征長度，具有較大的克努森數(shù)

Kn

，可以認為處在自由分子流的流動狀態(tài)。自由分子流狀態(tài)下的質(zhì)量流量為式(8):

式中，

n

為多孔板單孔個數(shù)，

?

為多孔板孔徑，m；

d

為水蒸氣流動長度，m；

m

為單個水分子質(zhì)量，kg；

k

為玻爾茲曼常數(shù)，

p

為真空環(huán)境壓力，Pa；飽和蒸氣壓

p

可根據(jù)式(9)給出:

此外，將式(8)帶入式(3)中，可得水升華器的升華熱流如式(10):

由于水不進入多孔板，此時水蒸氣的流動長度

d

即多孔板厚度

d

。聯(lián)立式(7)～(9)可以數(shù)值求解升華溫度

t

，并根據(jù)式(4)～(6)計算冰層厚度

d

、水層厚度

d

及底面溫度

t

。若熱載荷大于一定值，使升華表面蒸氣壓

p

大于水的三相點壓力

p

、升華溫度

t

超過三相點溫度

t

，水會進入多孔板轉(zhuǎn)換為周期模式。存在一個水升華器以升華模式工作的臨界熱載荷，標(biāo)志升華模式向周期模式的轉(zhuǎn)換，基于式(7)～(8)可得式(11):

當(dāng)水升華器的熱載荷小于式(11)時，認為處于升華模式；而當(dāng)熱載荷大于該值時，則認為水升華器處于周期模式下。此外，升華模式下熱載荷還需要滿足式(12)的限制，避免給水腔被冰層完全填充可能造成的損壞。這個最低熱載荷的數(shù)值同樣可以通過數(shù)值方法計算得出。

式中，

d

｜代表熱載荷

Q

下根據(jù)導(dǎo)熱關(guān)系式及式(6)計算得出的冰層厚度。

綜上，對于尺寸確定的水升華器，其升華模式穩(wěn)定工作的最高與最低熱載荷均是一個常數(shù)，可以通過式(11)～(12)確定。取升華面積0.04 m、多孔板厚度1 mm、給水腔厚度3 mm，可得圖5所示的結(jié)果，多孔板孔度、孔徑的增大會使熱載荷范圍變大、總體數(shù)值升高。

圖5 升華模式工作范圍Fig.5 Working scope of sublimation mode

2.3 熱載荷公式量綱分析

將式(11)中的小孔個數(shù)

n

用孔度

ε

展開并移項，兩邊再同除質(zhì)量流量

m

，可得式(13):

將上式進行量綱分析可以得到式(14):

式中，將式(13)的左邊將臨界熱載荷

Q

替換為實際熱載荷

Q

，可定義斯坦頓數(shù)

St

、雅各布數(shù)

Ja

倒數(shù)差為式(15):

式中，

ρ

為水蒸氣密度。定義為一個新的無量綱數(shù)

Su

如式(17)，包含水升華器的孔度、孔徑、面積等所有設(shè)計參數(shù):

結(jié)合密度、速度、面積比項，定義無量綱質(zhì)量流量如式(18):

根據(jù)2.2節(jié)的內(nèi)容，當(dāng)熱載荷大于一定值時水升華器將轉(zhuǎn)入周期模式。由于此前的分析保留了模式轉(zhuǎn)換對應(yīng)的三相點參數(shù)，因此熱載荷

Q

取臨界熱載荷

Q

時，可將式(11)轉(zhuǎn)化為如式(14)所示的無量綱形式，反映臨界熱載荷與水升華器設(shè)計參數(shù)之間的關(guān)系。對于孔度、孔徑及熱載荷分別在0.3～0.7、1～9μm、100～300 W間隨機取值，可得到如圖6所示的水升華器工作模式分布及模式轉(zhuǎn)換的無量綱特性曲線。對于水升華器不同設(shè)計參數(shù)組合的無量綱數(shù)

Su

，若熱載荷對應(yīng)的

St

取值在圖中的特性曲線以上，可以認為水升華器處于周期模式下工作，反之則認為進入了升華模式。此外，由式(12)決定的升華模式下邊界曲線同樣可以在圖中給出。

圖6 模式轉(zhuǎn)換無量綱特性曲線Fig.6 Dimensionless curve of mode transition

3 仿真模型

為了進一步探究水升華器升華模式的具體特性，基于2.2～2.3節(jié)理論分析，進行水升華器升華模式的數(shù)值建模與仿真計算。

3.1 焓-多孔模型

本文采用焓-多孔模型求解固-液相變問題。該模型以溫度和焓作為共同變量，在全計算區(qū)域內(nèi)建立統(tǒng)一形式的能量守恒方程。溫度和焓作為共同因變量的關(guān)系可以式(19)用來表示:

引入糊狀區(qū)(Mushy zone)的概念，將溫度處于固相線溫度與液相線溫度之間的區(qū)域視為糊狀區(qū)，該區(qū)域內(nèi)網(wǎng)格的液相分?jǐn)?shù)視為孔度，通過不斷更新所有網(wǎng)格單元內(nèi)的孔度追蹤固液界面的位置。液體分?jǐn)?shù)的表達為式(20):

對于水-冰相變問題，其固相線與液相線溫度相等，均為273.15 K。糊狀區(qū)的動量表達為式(21):

式中，

β

為液體分?jǐn)?shù)，

ε

為防止分母為0的小量，

A

為糊狀區(qū)常數(shù)，

V

為速度矢量，

V

為固態(tài)離開原區(qū)域的漂移速度?；谑?20)～(21)，可得相變問題的控制方程組如式(22)～(24)所示。

連續(xù)性方程為:

由于水升華器工作在失重的空間環(huán)境中，可以忽略重力的影響，得到動量方程:

能量方程為:

式中，總焓

H

為顯焓

h

與潛熱焓Δ

H

之和，關(guān)系如式(25):

式中，

h

為參考焓值，J/kg；

T

為273.15 K的參考溫度；

L

為水-冰相變潛熱，J/kg。

3.2 網(wǎng)格劃分及Fluent求解模型

使用商業(yè)軟件ANSYS ICEM CFD及ANSYS Fluent 2020 R1版本作為網(wǎng)格劃分工具及求解器。由于供水流量很小，流動模型選擇層流Laminar。能量方程的求解調(diào)用Fluent中基于焓-多孔模型的Solidification/Melting模型，選擇基于壓力的隱式穩(wěn)態(tài)單精度求解器進行求解，當(dāng)整體熱量達到平衡后計算結(jié)束。

仿真分析相同外形的水升華器在不同孔度、孔徑時的升華模式特性，水升華器結(jié)構(gòu)參數(shù)及相關(guān)常數(shù)如表1所示，其中3個工況已在圖7中標(biāo)出。

表1 水升華器參數(shù)及相關(guān)常數(shù)Table 1 Parameters of sublimator and relative constants

對給水腔進行二維矩形結(jié)構(gòu)網(wǎng)格的劃分，并且選取了數(shù)量為7150、8400、10500的網(wǎng)格進行仿真網(wǎng)格無關(guān)性驗證計算。3套網(wǎng)格的工況2計算結(jié)果如圖7所示，不同網(wǎng)格數(shù)量下的計算結(jié)果最大差異不超過0.15%。

圖7 網(wǎng)格無關(guān)性驗證Fig.7 Grid independence validation

綜合考慮計算的準(zhǔn)確性及計算資源的占用，最終選取600×14共8400數(shù)量的網(wǎng)格，網(wǎng)格劃分及邊界條件設(shè)置情況如圖8所示。其中，底面熱載荷所在的加熱面設(shè)置為等熱流密度邊界，為根據(jù)圖6確定的熱載荷區(qū)間選取的合適的熱載荷。頂面為升華界面，定義為由溫度決定的熱流邊界，與孔度、孔徑及升華溫度等有關(guān)，如式(9)～(10)所示。計算域右側(cè)加入一個由式(7)確定的等質(zhì)量流量入口，同時在頂面施加一個質(zhì)量源項UDF模擬升華帶走的質(zhì)量流。邊界條件結(jié)果如表2所示。

圖8 網(wǎng)格與邊界條件Fig.8 Mesh and boundary condition

表2 不同工況邊界條件取值Table 2 Boundary condition settings

4 結(jié)果與分析

由于Fluent中是基于式(26)計算液體分?jǐn)?shù)的，松弛因子

γ

的存在會導(dǎo)致計算域中出現(xiàn)糊狀區(qū)。為避免不確定性，采用

t

=273.15 K等溫面表示水-冰界面如式(26):

三工況整體分布較為相似，取工況2的計算結(jié)果進行分析。整體及入口附近液體分?jǐn)?shù)云圖如圖9所示，紅色代表液態(tài)水，綠色、藍色代表相變形成的糊狀區(qū)及固態(tài)冰?？梢娛艿綔囟容^高的供水的影響，入口附近幾乎沒有形成冰層；而向后直至給水腔盡頭，入口效應(yīng)逐漸消失，形成了厚度較為一致且水-冰界面與上下邊界大致平行的薄冰層。

圖9 液體分?jǐn)?shù)云圖Fig.9 Contour of liquid fraction

工況2的流場矢量圖如圖10所示。受固液界面位置的影響，水流從入口進入后僅在下側(cè)的液態(tài)區(qū)域中流動，上側(cè)的固態(tài)冰層區(qū)域中不存在流動；隨后水流繼續(xù)向前流動遠離入口，速度逐漸降低并最終停止流動。此外，還能觀察到在冰面及底部壁面附近的流動邊界層效應(yīng)。

圖10 速度矢量圖Fig.10 Vector of velocity

工況2的溫度分布云圖如11所示。可以發(fā)現(xiàn)，溫度分布與液體分?jǐn)?shù)分布具有很強的關(guān)聯(lián)性。同樣將入口附近放大，可見入口的高溫水流前進的同時很快地受到了冷卻，等溫線逐漸傾斜；入口效應(yīng)逐漸消失后，受到上下熱流邊界的影響，等溫線也逐漸平行于上下邊界。此外，還能觀察到右下角區(qū)域在底面加熱作用下達到了全域最高溫度，已超過了293.15 K的供水溫度。

圖11 溫度云圖Fig.11 Contour of temperature

統(tǒng)計各工況計算穩(wěn)定后的底面溫度、升華溫度與其對應(yīng)的冰層厚度的仿真結(jié)果與理論計算值對比，如表3所示。各工況的溫度結(jié)果相對誤差都在0.5%以內(nèi)，底面溫度仿真結(jié)果整體比理論值高1 K左右，而升華溫度仿真結(jié)果比理論值稍低。造成這種誤差的原因是計算域頂部升華熱流、升華溫度耦合的邊界條件設(shè)置，其熱流是由升華界面平均溫度決定的，因此入口高溫水流對升華界面的邊界有著不可忽略的影響。

表3 溫度、界面位置結(jié)果對比Table 3 Comparison of temperature and thickness of ice

各工況冰層厚度的計算結(jié)果與理論值相對誤差的絕對值均小于1%。同樣受到入口效應(yīng)的影響，計算域右側(cè)沒有形成冰層，同時溫度的誤差也會影響冰層厚度，因此冰層平均厚度整體較理論值稍低，對應(yīng)的相對誤差偏向于負值。

5 結(jié)論

1)分析推導(dǎo)得到了水升華器工作模式轉(zhuǎn)換的特性曲線，該曲線能給出不同熱載荷及孔度、孔徑等結(jié)構(gòu)參數(shù)下水升華器的工作模式；

2)對不同孔度、孔徑共3種工況進行仿真并與理論值比較，各溫度結(jié)果的相對誤差均在0.5%以內(nèi)、冰層厚度的相對誤差絕對值均在1%以內(nèi)，仿真模型通過了理論驗證；

3)入口附近較高溫度的供水會造成等溫線傾斜、無冰層形成等入口效應(yīng)及冰層厚度仿真結(jié)果偏低等誤差。

本文僅研究了低熱載荷下的升華模式，可在后續(xù)研究中開展對熱載荷較大時的周期模式的理論分析及數(shù)值仿真工作。