基于相場法熔鹽相變的模擬研究

王超穎，顧璠

（東南大學能源與環(huán)境學院，江蘇南京 210018）

0 引言

利用相變蓄能技術可對不連續(xù)、不穩(wěn)定的熱量進行充分利用，有效提高能源利用效率[1]，在太陽能熱利用[2]、工業(yè)余熱回收[3]、建筑節(jié)能等領域中相變材料都有著誘人的應用前景。陳思等[4]比較了不同的蓄熱熱泵，得出制冷劑化霜時先流過熱交換器再流過蓄熱器更合理的結論；王聰等[5]研制了一種新型的管肋式相變蓄冷換熱器，經測試蓄冷器的 平均功率為4,057.2 W，為太陽能空調實驗系統(tǒng)提供了設計依據。KALAISELVAM等[6]研究了相變材料在柱狀密封體中完全相變所用的時間，通過模擬和實驗對比可知完全凝固時間依賴于Stefan數和生成熱參數，完全熔化時間依賴于當量導熱系數，模擬熔化過程比凝固過程更接近實驗結果。由于熔鹽材料在蓄熱技術中的廣泛應用，本文擬用相場法模擬計算其相變過程。

相場法基于Ginzburg-Landau相變理論[7]，采用擴散界面模型，認為兩相界面有一定厚度，各種傳輸現象均通過該界面區(qū)域完成。以熔化過程為例，固相、液相和界面作為不同的相態(tài)用體積分數來跟蹤系統(tǒng)中的相，無需追蹤界面將固液相分開處理，因此在整個求解域中均采用相同的數值計算方法。相場法早期用于枝晶生長的研究[8]，BECKERMANN等[9]用一種相場模型模擬了凝固過程中的微觀結構變化，如今，相場法已成為許多領域一種重要的研究方法[10-11]。

1 計算方法與驗證

1.1 理論模型

相場法用相場變量Φ表示系統(tǒng)時間和空間上的物理狀態(tài)。在本文的計算中，Φ的值代表了液態(tài)或固態(tài)，Φ=1表示固態(tài)，Φ=-1表示液態(tài)，-1＜Φ＜1代表相界面。

流體的速度場和壓力場通過不可壓N-S方程表示[12-13]：

對相場變量的處理方程[14]：

通過相界面的質量(相變速率)計算如下[15]：

連續(xù)性方程考慮相變后變更如下：

式中：F為體系自由能；γ為遷移率；λ為混合能量密度；ΔHl為熔化焓。

1.2 相方程用于液固相變計算的討論

假定初始物理場溫度一致，處于相變溫度，考慮重力。流動中，認為邊界無滑移。圖1為二維等壓流動過程入口流體速度積分，可知當黏性系數大于1 Pa·s時，速度極緩且趨于穩(wěn)定，在之后的相變計算中，固相的黏度取1 Pa·s。

圖1 二維等壓流動過程入口流體速度積分

1.3 驗證計算

本文查閱已有的實驗結果[16]，模擬計算來驗證計算模型的可靠性。石蠟物性如下：石蠟的熔點為41.98 ℃，相變潛熱為141 kJ/kg，導熱系數為0.2 W/(m·K)，比熱為2.0 kJ/(kg·K)，液態(tài)密度為800 kg/m3，固態(tài)密度為916 kg/m3，液體運動黏度為5×10-6m2/s。

驗證計算模擬了石蠟在方形蓄熱槽中的相變過程，蓄熱槽中心有一圓管，管內提供60 ℃的恒溫水，石蠟分布在管外，水與石蠟通過管壁完成換熱，方形槽邊長為10 cm，圓管直徑為2 cm。

圖2為石蠟測點分布，上測點1#到管中心距離為1.5 cm，下測點2#到管中心距離為4.5 cm。圖3為蓄熱槽的結構。為了簡化計算模型假設：管壁材料的導熱系數很大，管壁厚度可視為零，管壁內側與外側溫度相等；相變材料在相變時不發(fā)生過冷和過熱現象；流動為二維，非穩(wěn)態(tài)，層流，蓄熱槽沒有漏熱損失，即蓄熱槽上下左右表面絕熱。

圖2 石蠟測點分布

圖3 蓄熱槽的結構

圖4為上測點和下測點溫度的對比。蓄熱槽內，上測點處的相變材料溫度上升迅速，在短時間內完成相變并達到穩(wěn)定狀態(tài)，且上部液態(tài)區(qū)形成了自然對流加快了傳熱速度；下測點處的相變材料溫度變化緩慢，在300 min內都沒有達到相變溫度，一是因為石蠟導熱系數較小，二是因為下部區(qū)域傳熱方式為導熱。從文獻中的實驗溫度值與本研究的模擬溫度值對比結果可看出，對初始溫度的假設并不影響最終結果，模擬與實驗達到穩(wěn)定相變的時間相近，最終溫度吻合良好，驗證了相場模型計算的可靠性。

圖4 上測點和下測點溫度的對比

2 二元系液固相變模擬和特性分析

由于石蠟蓄熱能力較小，選取混合鹽作進一步的研究。混合鹽經適當改變其組分的配比可得到希望的熔點，溫度適用范圍更廣，將儲熱性能好的高價格物質與低價格物質合在一起使用可節(jié)省成本[17]。綜合經濟因素與導熱性能，本文選用了NaCl-CaCl2作為相變材料，研究高溫熔鹽的相變特點。圖5為某實驗的NaCl-CaCl2二元相圖[18]，當NaCl的摩爾占比約為50%時，體系的共晶溫度為773 K。

混合工質物性按摩爾比計算，工質物性參數如下[19-20]：二元系的共熔點為773 K，相變潛熱為377.303 kJ/kg，導熱系數為2.61 W/(m?K)，比熱為1.07 kJ/(kg·K)，液態(tài)密度值為1,861 kg/m3，固態(tài)密度值為2,159.6 kg/m3，液體運動黏度為1.2×10-6m2/s。利用上一節(jié)介紹的模型，初始溫度為相變溫度，模擬材料在不同結構、不同壁溫下的熔化過程。

圖5 氯化鈉與氯化鈣二元相圖

2.1 液固相空間分布與傳熱特性

對蓄熱槽內熔鹽的熔化過程進行模擬，管內提供1,200 K的恒溫熱源，分別得到管壁熱流密度和相界面長度的變化曲線（圖6～圖7）。

圖6 管壁熱流密度變化曲線

圖7 相界面長度變化曲線

熱量從管壁往四周輻射，隨著周圍溫度上升，溫差變小，且隨著相變的進行，液相增多，相界面長度（固液接觸長度）不斷變長，至26 s時達到最大，同時相變熱大量積聚。管壁熱流量的值開始時特別高，原因是管壁與很薄的液體層接觸，溫差很大，熱阻很?。浑S著相變進行，管壁周圍液體層逐漸變厚，傳熱熱阻變大，熱流量持續(xù)減?。辉?0 s附近熱流量增大，是因為上部液體層出現自然對流，加快了熱量傳遞，此時相變區(qū)域覆蓋到將近整個區(qū)域的一半。而后由于固相變少，固液接觸減少，界面長度不斷縮短。當上部區(qū)域基本熔化、傳熱向下方擴散時，相界面發(fā)生變化，下方的低溫區(qū)傳熱方式為導熱，熔化緩慢，熱流密度保持降低趨勢直到熔化完成。

2.2 不同壁溫對比

圖8為圓管內測點分布，外管徑為4 cm，內管徑為1 cm，上部測點距圓管中心0.8 cm，下部測點距圓管中心1 cm。圖9為上測點和下測點溫度變化圖，由圖9可知，當壁溫Tw=1,400 K時，圓管上部分相變材料溫度迅速升高，在短時間內達到一個相對穩(wěn)定的傳熱過程，壁溫1,200 K相對1,400 K溫度升高慢，但趨勢類似，都是在短時間內有一個較大的升溫，Tw=1,000 K時，斜率比較低，溫度上升比較慢。而圓管下部分熔化速度緩慢，經過60 s，1,400 K壁溫下的該測點材料才開始熔化，Tw=1,200 K時用了90 s，而Tw=1,000 K時，在120 s之前幾乎沒有開始熔化。綜上所述，上部區(qū)域形成對流，傳熱速度快，下部傳熱方式為導熱，傳熱速度慢。因此，在利用相變蓄熱的結構中，應注意強化下方區(qū)域的傳熱性能。

2.3 不同結構對比

圖10所示為當Tw=1,200 K時，圓管內相分數的變化。由圖10可知，開始時的熔化過程主要集中在圓管上部區(qū)域。上部區(qū)域內的相變材料快速熔化，形成對流，當液相局部布滿頂端時，熔化過程向左右兩端擴散，形成類似折扇的圖案，而后向下方繼續(xù)擴展開來。但側面熔化完成后，熔化速度減慢，此時的傳熱方式為導熱，在下方區(qū)域的熔化過程（時間幾乎為上方熔化時間的兩倍）中，由于表面張力的影響，相變材料沿著管壁在圓管底部形成凸起，最終慢慢熔盡。

圖8 圓管內測點分布

圖9 上測點和下測點溫度的變化

圖10 圓管內相分數的變化（Tw=1,200 K）

圖11為相同壁溫下方管內的相分數變化，方管與圓管的熔化過程略有不同，方管耗時時間短。一方面，外管周長相同，圓管相對于方管傳熱面積大，所需的熔化時間更長。另一方面，由于形狀不同，方管內的熔化過程受到的表面張力影響比較小，雖然熔化至相界面將近水平時耗時相近（約40 s），但熔化后期方管內的傳熱面積基本保持不變，傳熱速率穩(wěn)定；而圓管內壁面效應積聚，底部固相區(qū)遲遲難以熔化，傳熱速率極大地減小。

圖11 方管內相分數的變化（Tw=1,200 K）

3 結論

本文模擬了石蠟的熔化過程并將其與文獻中的實驗結果對比，驗證了相場模型用于液固兩相流相變計算的可靠性。研究了二元相變材料NaCl-CaCl2在圓形與方形蓄熱槽內伴有對流傳熱的蓄熱過程，給出了定壁溫下熔化過程的固液相圖，得到如下結論：

1）壁溫越高，熔化越快。相變過程中的傳熱性能受黏度的影響，黏度增大（黏度大于1 Pa?s時基本穩(wěn)定）使相界面更平穩(wěn)但相應地傳熱速度會減慢，因此要根據對蓄熱過程的要求選擇合適的材料；

2）從相界面長度的變化及傳熱特性的分析可知，熔化前期用時短，因為對流引起傳熱系數增大，而熔化后期速度慢，傳熱方式為導熱；

3）相同周長的圓管比方管耗時多，且方管受到表面張力的影響比圓管小。