郭茶秀，陳 俊

（鄭州大學(xué) 化工與能源學(xué)院，鄭州 450001）

方腔內(nèi)合金相變材料的熔化過程*

郭茶秀?，陳 俊

（鄭州大學(xué) 化工與能源學(xué)院，鄭州 450001）

本文利用CFD軟件對方腔內(nèi)合金相變材料在溫差作用下的熔化過程進行了數(shù)值模擬研究。通過與試驗結(jié)果的比較，驗證了本文采用的模型和算法的正確性。同時，詳細(xì)研究了熔化過程中合金相變材料的固?液界面、溫度及環(huán)流速度的變化規(guī)律。結(jié)果表明，為了減小傳熱熱阻，加快合金的熔化速度，從方腔的側(cè)面和底面加熱相變材料時熔化效果最好。

合金相變材料；熔化；數(shù)值模擬

0 引 言

工貿(mào)企業(yè)存在大量的低品位工業(yè)余熱，由于工業(yè)余熱的不連續(xù)性，使得現(xiàn)有以直接回收利用為主的余熱回收技術(shù)無法實現(xiàn)低溫的余熱跨時間、空間梯級利用，從而導(dǎo)致企業(yè)能源綜合利用率較低。而儲能技術(shù)能解決能量供求在時間和空間上分配的不均勻性，可有效提高能源利用率，在工業(yè)的熱過程、新能源利用、民用建筑采暖等方面已得到廣泛應(yīng)用[1]。

目前的熱能存儲技術(shù)主要包括顯熱、潛熱和化學(xué)儲熱三種。其中潛熱存儲以其儲熱密度高、溫度變化小、容易控制等特點，成為合理利用能源、減小環(huán)境污染、提高能源利用率的一種儲能方式[2]。潛熱存儲所使用的相變材料分為有機、無機和混合相變材料三類，其中Ga、Pb等合金，因具有儲能密度大、熱循環(huán)穩(wěn)定性好、導(dǎo)熱系數(shù)大、過冷度小等特點而受到國內(nèi)外研究者的關(guān)注，并且主要在合金的熱物理性能研究方面取得了一定的進展[3]， 但在其傳熱性能方面研究較少，因此在合金的傳熱方面還可以開展許多研究。

在合金相變材料相變過程中，由溫差引起的自然對流始終存在，使得合金相變材料熔化或凝固時固?液界面的位置受時間、溫度、傳熱量以及容器形狀尺寸等因素的影響，極大地增加了研究固?液界面運動規(guī)律和傳熱效率的復(fù)雜性。本文利用CFD軟件FLUENT[4]對合金相變材料在自然對流條件下的相變過程進行了數(shù)值模擬研究，探討了其相變過程中自然對流對溫度分布、環(huán)流速度及固?液界面隨時間的影響，為相變材料的有效利用、蓄熱系統(tǒng)的優(yōu)化設(shè)計及性能的準(zhǔn)確分析提供參考依據(jù)。

1 物理模型及數(shù)學(xué)模型的建立

1.1 物理模型

本文研究的物理模型[5]如圖1所示，方腔內(nèi)裝入合金鎵作為相變材料，其熱物性如表1所示。方腔的尺寸為：x方向長8.89 cm，y方向長6.35 cm，z方向遠(yuǎn)大于x、y方向。上下兩面設(shè)有保溫層，左右兩側(cè)面溫度恒定，相變材料的相變溫度在兩者之間以便于模擬計算。左側(cè)面為加熱面，溫度為38℃；右側(cè)面溫度為28.3℃。為了分析方便，可將其物理模型簡化為x、y平面內(nèi)的二維模型，并做了以下假設(shè)[6]：

（1）流體中粘性耗散不計；（2）除密度外其他物性為常數(shù)；（3）自然對流是由液相的密度發(fā)生變化而形成的，而密度變化符合Boussinesq假設(shè)。

圖1 物理模型Fig. 1 Physical model

表1 相變材料熱物性參數(shù)Table 1 property parameters of PCM

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

由于合金相變材料相變過程是一個非穩(wěn)態(tài)過程，且相變過程中傳熱方式以導(dǎo)熱和對流為主，比單純的導(dǎo)熱或者單純的對流情況復(fù)雜得多。尤其是因為相變傳熱是一個非線性的傳熱問題，疊加原理不再適合。因此在求解能量方程的過程中引入焓，并將整個相變區(qū)域看成多孔介質(zhì)，孔隙度用液相分?jǐn)?shù)lf來表示，當(dāng)lf=0時，模擬對象為固相，流體的流速為0；當(dāng)lf=1時，模擬對象為液相，流體的流速為動量方程的速度；當(dāng)0＜lf＜1時，模擬對象為模糊區(qū)，流速為液相分?jǐn)?shù)的函數(shù)?？紤]到自然對流的影響，在連續(xù)方程中應(yīng)考慮到密度的變化，在動量方程中應(yīng)該考慮體積力的影響，故其二維控制方程及邊界條件如下：

連續(xù)方程：

動量方程：

能量方程：

基于前面的假設(shè)，動量方程中的密度表示為：

動量方程的源項：

能量方程源項：

式（6）中ε是一小于0.000 1的數(shù)，Amush為糊狀區(qū)域的連續(xù)數(shù)（在104～107之間）。

能量基本方程式（4）中的H表示為：

自然對流是在重力場或者其他力場的作用下由密度差引起的浮升力產(chǎn)生的，大多數(shù)情況下，密度差是由于溫度差引起的，液體受熱或冷卻后體積膨脹或縮小，導(dǎo)致密度增大或減小。假設(shè)遠(yuǎn)離加熱面某點處的液體密度為ρ'，加熱面附近某點處液體的密度為ρ，則ρ=ρ'(1-βΔT)，于是這兩點將形成壓差：

當(dāng)Δt較小時，，在壓差的推動下，必然造成液體的環(huán)流，其速度滿足關(guān)系式：

式中，β為體積膨脹系數(shù)，ΔT為兩點的溫差，l為受熱面的高度，g為重力加速度，u為速度絕對值，ρ為密度，p為壓力。

可見環(huán)流速度的絕對值取決于流體的性質(zhì)（例如粘度）、流動空間的幾何形狀及尺寸。

2 模擬過程的實現(xiàn)及結(jié)果

2.1 參數(shù)設(shè)置

利用Gambit 2.1.3對圖1所示簡化后的2D計算模型進行網(wǎng)格劃分，經(jīng)過網(wǎng)格獨立性考核，網(wǎng)格總數(shù)選為1 344個。將網(wǎng)格文件導(dǎo)入FLUENT 6.3中，定義求解器Solver為分離式2D，非穩(wěn)態(tài)，一階隱式；選定Solidification/Melting模型；輸入相變材料熱物性參數(shù)及邊界條件；離散項中壓力梯度項采用PRESTO格式，壓力速度耦合項采用SIMPLE算法；時間步長取0.1 s。為了得到相變過程中溫度與速度的關(guān)系，在計算模型中設(shè)置了水平線y=3.175上的4點：point1、point2、point3、point4。4個點的橫坐標(biāo)為x=1.778 cm、x=3.556 cm、x=5.334 cm、x=7.112 cm，并監(jiān)視這4個點的溫度、速度與時間的關(guān)系。另外，為了獲得各個時刻系統(tǒng)中溫度的分布情況，定義了平行于x軸的4條直線，這4條直線的縱坐標(biāo)為y=1.27 cm、y=2.54 cm、y=3.81 cm、y=5.08 cm[5]。

2.2 模擬結(jié)果及分析

圖2所示為模擬計算結(jié)果，與文獻(xiàn)[7]的試驗結(jié)果相比較，證明本文模擬的結(jié)果是正確的。

圖2 計算與實驗結(jié)果的比較Fig. 2 Comparison of simulation and experimental results

為了分析蓄熱溫差ΔT=Thot-T（cTc為相變材料的熔點）為10℃時自然對流對相變過程的影響規(guī)律，將初始溫度設(shè)為相變材料的熔點。圖3所示為ΔT=10℃，計算模型中相變材料的固?液界面在不同時刻的固?液界面分布情況。

圖3 相變材料固?液界面Fig. 3 Interface of solid-liquid in PCM

從圖3可以看出，在熔化的初始階段，由于相變材料的液態(tài)量較少，傳熱的方式主要以導(dǎo)熱為主，2 min以前方腔內(nèi)的界面曲線彎曲不大，但是隨著時間的推移，固?液界面彎曲越來越嚴(yán)重。由于液態(tài)量逐漸增多，導(dǎo)熱方式從以導(dǎo)熱為主轉(zhuǎn)變成以對流為主，靠近加熱壁面的液相相變材料由于溫度升高，密度減小，在浮力的作用下向上移動，方腔的上部受到熱流的沖擊，使得上半部分相變材料熔化的速度加快，同時由于熱的相變材料在與固?液界面發(fā)生對流換熱后，液態(tài)的溫度降低，使得液態(tài)相變材料在靠近固?液界面區(qū)域向下流動，產(chǎn)生回流。由于沿固?液界面回流時液態(tài)溫度逐漸降低，故使得固?液界面的移動速度沿自身自上而下逐漸減小。

圖4 和圖5為自然對流條件下，相變材料分別在5 min、15 min時刻的相變過程中的等溫線圖。從不同時刻等溫線來看，受自然對流的影響，在液相區(qū)（溫度大于302.93 K）有兩個區(qū)域的溫度梯度較大，一個是位于方腔上部的固?液界面前沿，此區(qū)域是由于受到熱流體的沖擊；另一個區(qū)域是位于方框下方靠近加熱壁的區(qū)域，此區(qū)域是由于冷流體回流后被加熱。受自然對流的影響，一方面，液體加熱后溫度升高，密度減小，在浮力的作用下向上流動，加熱面附近的液態(tài)相變材料的溫度迅速升高并沿加熱面上升；另一方面，在固?液界面附近的相變材料，由于不斷地向固?液界面?zhèn)鬟f熱量，溫度不斷降低并導(dǎo)致密度增大，逐漸下沉。高溫流體在上升過程中與固?液界面接觸時，由于兩者的溫差較大（固?液界面的溫度始終等于相變溫度），使得該區(qū)域的溫度梯度較大；同樣低溫流體下降冷卻后與加熱面接觸時，由于加熱面的溫度遠(yuǎn)高于流體的溫度，使得該區(qū)域的溫度梯度也較大。低溫流體加熱后上升，高溫流體傳遞熱量后下降，從而導(dǎo)致液態(tài)的循環(huán)流動。因此，可以推測從側(cè)面和底面加熱時，固態(tài)能夠始終與加熱面相接觸，熔化速度最快，效果最好。

圖4 相變材料在5 min時的等溫線Fig. 4 Isothermal line at 5 minutes in PCM

圖5 相變材料在15 min時刻的等溫線Fig. 5 Isothermal line at 15 minutes in PCM

圖6為=ΔT10℃條件下，本文所設(shè)置的4點溫度隨時間的變化圖。由圖6可看出，未達(dá)相變溫度之前，各點溫度逐漸升高，距離加熱面越近，曲線斜率越大，溫度升高越快；當(dāng)溫度達(dá)到相變溫度之后，各點的溫度發(fā)生波動，波動先增大后減小，且溫度始終不低于相變溫度。溫度未達(dá)相變溫度之前，相變材料為固態(tài)，熱量以導(dǎo)熱的形式傳遞，離加熱面越近，熱阻越小，溫度上升越快，此階段以顯熱的形式儲存熱量；當(dāng)溫度到達(dá)相變溫度之后，受液態(tài)自然對流的影響，該點溫度產(chǎn)生波動，此階段以潛熱的形式儲存熱量；當(dāng)相變完成后，自然對流引起的循環(huán)流動使流體區(qū)域的溫度均勻化，隨著蓄熱的進行，溫度波動逐漸減小并逐漸升高，此階段又以顯熱的形式儲存熱量。在不考慮自然對流的條件下，各點的溫度均不會發(fā)生波動[8]。

圖6 四點溫度隨時間的變化Fig. 6 Temperature changes with time at four points

圖7～圖10分別為y=1.27、2.54、3.81、5.08 cm四條直線上，合金相變材料在2、4、6、8、10、12、14、16 min時的溫度分布圖。從圖中可以看出在2 min的時候，相變材料在相變的初始階段，由于液相量比較少，各等高線上的溫度在熔點以上呈直線分布；4 min的時候，由于熱流體向上流動，y=1.27 cm和y=2.54 cm 的溫度仍呈線性分布，而y=3.81 cm和y=5.08 cm的溫度分布線發(fā)生了扭曲，并隨時間逐漸加強；16 min的時候，在x=0.01 m到x=0.08 m的中間部分，直線y=5.08 cm的溫度分布比較平穩(wěn)、溫度最高，y=3.82 cm和y=2.54 cm次之，y=1.27 cm最陡。蓄熱過程中，液態(tài)相變材料加熱后不斷上升，并盡量靠近上壁面流動，溫度越高，上升越高，使同一水平面上的溫度均勻化。由以上可知：在無限接近上底面的等高線上，溫度分布最平穩(wěn)，最接近加熱面的溫度。合金相變材料在相變過程中受自然對流的影響，液態(tài)區(qū)域的傳熱以對流為主，液態(tài)受熱之后密度減小，在浮升力的作用下迅速離開原來的位置，使得熱量盡量傳遞到溫度較低區(qū)域，導(dǎo)致固態(tài)熔化速度加快。

圖7 y=1.27 cm溫度隨時間分布Fig. 7 Temperature distributions at liney=1.27 cm

圖8 y=2.54 cm溫度隨時間分布Fig. 8 Temperature distributions at liney=2.54 cm

圖9y=3.81 cm溫度隨時間分布Fig. 9 Temperature distributions at liney=3.81 cm

圖11所示分別為point1、point2、point3、point4四點的速度隨時間的變化規(guī)律。結(jié)合圖7可以看出：相變材料在未達(dá)熔點之前，各點速度為零。當(dāng)溫度超過熔點以后，由于對流產(chǎn)生速度的影響，各點的溫度在熔點以上發(fā)生波動。同時由于環(huán)流的流速與溫差成正比，各點的速度隨溫度波動增大而增大，到達(dá)最大之后逐漸減小。隨著溫差波動的減小，各點速度逐漸趨于零。

圖10 y=5.08 cm溫度隨時間分布Fig. 10 Temperature distributions at liney=5.08 cm

圖11 四點速度隨時間的變化Fig. 11 Velocity changes with time at four points

圖12為溫差ΔT=10℃、15℃，相變時間為4 min時方腔內(nèi)的相變材料的速度矢量分布圖。從圖中可以看出：在同一時刻，溫差為15℃時相變材料的最大速率為2.36 cm/s，而溫差為10℃時的最大速率為1.74 cm/s，即溫差越大，相變材料的液相速率越大。在液相區(qū)，受自然對流的影響，液相環(huán)流隨液相所占比例的增大而越來越明顯，環(huán)流速度的絕對值除了取決于冷熱面的溫差外，還與流體的性質(zhì)（如粘度、膨脹系數(shù)等）、流動空間的幾何形狀尺寸有關(guān)。

圖12 速度矢量分布圖（4 min時刻）Fig. 12 Velocity vector distributions at 4 minutes

3 結(jié) 論

（1）本文對自然對流作用下的合金相變材料鎵的熔化過程建立了數(shù)學(xué)模型和物理模型，并利用FLUENT軟件對其進行了數(shù)值模擬。通過與試驗結(jié)果的比較，表明本文對合金相變材料的相變問題進行傳熱分析是正確的。

（2）合金相變材料在熔化過程中，由于浮升力的作用，熱流體上升，冷流體下降，沿固?液界面回流時液態(tài)溫度逐漸降低，故使其固?液界面的移動速度沿自身自上而下逐漸減小。在自然對流作用下，環(huán)流速度與溫度差成正比。環(huán)流速度隨溫度波動增大而增大，到達(dá)最大之后逐漸減小，隨著溫差波動的減小，速度逐漸趨于零[5]。

（3）相變材料水平線的位置越高，溫度分布越平穩(wěn)，由此可推出無限接近上表面的水平線上，溫度分布最平穩(wěn)，最接近加熱面的溫度。

（4）根據(jù)對自然對流對合金相變材料熔化過程影響的分析，可知從方腔的兩個側(cè)面和底面加熱相變材料時熔化效果最好，因為在相變材料熔化過程中固態(tài)下沉且更接近加熱面，從而減小傳熱熱阻，使熔化速度加快。

[1] Agyenim F, Hewitt N, Eames P, et al. A review of materials, heat transfer and phase change problem formulation for latent heat thermal energy storage systems, Renewable and sustainable energy reviews[J]. 2010, 14(2): 615-628.

[2] 張寅平, 胡漢平, 孔祥東. 相變貯能——理論和應(yīng)用[M]. 中國科技大學(xué)出版社, 1996.

[3] 李元元, 程曉敏, 俞鐵銘. Sn-Zn合金相變儲熱材料的熱循環(huán)穩(wěn)定性[J]. 特種鑄造及有色金屬, 2012, 32(7): 674-676.

[4] 韓占忠, 王敬, 蘭小平. FLUENT流體工程仿真計算事例與應(yīng)用[M]. 北京理工大學(xué)出版社, 2004.

[5] 陳俊. 相變蓄熱技術(shù)的數(shù)值模擬研究[D]. 鄭州: 鄭州大學(xué), 2007.

[6] 陶文銓. 數(shù)值傳熱學(xué)[M]. 西安: 西安交通大學(xué)出版社, 2001.

[7] Gau C, Viskanta R. Melting and Solidification of a Pure Metal on a Vertical Wall[J]. ASME Journal of Heat Transfer, 1986, 108(1): 174-181.

[8] 周建輝. 高溫固液相變蓄熱器的結(jié)構(gòu)設(shè)計和數(shù)值模擬[D]. 哈爾濱: 哈爾濱工程大學(xué). 2004.

Phase Change Process of Alloy Phase Change Material in a Rectangular

GUO Cha-xiu, CHEN Jun
(School of Chemical and Energy Engineering, Zhengzhou University, Zhengzhou 450001, China)

The melting process of alloy phase change material in a rectangular under temperature difference is simulated by CFD software. The simulation model is verified by the experimental result. The impacts of the natural convection on the solid-liquid interface, the temperature and the circulation speed in the melting process are analyzed. The results show that the phase change material in rectangular should to be heated from the two sides and the bottom in order to reduce heat transfer resistance and get a faster melting speed.

alloy phase change material; melting; numerical simulation

TK02

A doi：10.3969/j.issn.2095-560X.2015.04.012

2095-560X（2015）04-0319-06

郭茶秀（1968-），女，博士，教授，主要從事新能源的開發(fā)與利用方面的研究。

2015-04-14

2015-05-28

國家自然科學(xué)基金（51176173）；河南省教育廳科學(xué)技術(shù)研究重點項目（14A480002）

? 通信作者：郭茶秀，E-mail：guochaxiu@163.com