基于二維梯度樹狀肋相變儲熱系統強化傳熱機理

張欣宇，楊曉宏，2，張燕楠，徐佳錕，郭梟，田瑞，3

（1 內蒙古工業大學能源與動力工程學院，內蒙 古呼和浩特 010051； 2 風能太陽能利用技術教育部重點實驗室，內蒙古呼和浩特 010051； 3 內蒙古可再生能源重點實驗室，內蒙 古呼和浩特 010051）

引 言

在能源應用中，太陽能由于清潔無污染、綠色環保、輻射能流高、可再生而在全球廣泛應用，但其應用受到季節、陰雨、晝夜更替等自然條件影響。為解決太陽能在時間、空間上不匹配的問題[1-2]，采用儲能技術將暫時不用或多余的太陽能以熱能的形式存儲起來[3-4]。相變儲熱在合適的溫度范圍內具有較高的能量存儲密度，且過程近乎恒溫，具有良好的熱穩定性，被廣泛應用[5-7]。但相變材料（phase change material,PCM）的導熱性能差[熱導率普遍低于0.5 W∕(m·K)]，直接影響了相變儲熱系統能量存儲密度和傳熱速度[8]。研究發現，強化相變儲熱系統可以改變相變儲熱系統的結構，通過肋片增大傳熱流體（heat transfer fluid,HTF）與PCM 的換熱面積[9-11]；同時也可以制備高熱導率的復合相變材料，通過向基礎相變材料中添加金屬粒子或高導熱性能的骨架[12-13]實現熱流體與相變材料的傳熱強化。

通過合理的肋片設計來提高相變儲熱系統溫度分布的均勻性是強化傳熱的有效途徑[14-15]。Khan等[16]研究了水平潛熱儲能單元中肋片方向對強化傳熱的作用。基于焓孔隙率方法，對具有縱向肋片的管殼式相變儲熱系統的二維瞬態模型進行數值計算，通過改變翅片的不同角度（0°≤θ≤90°）進行了研究。Yu 等[17]設計了梯度肋片相變儲熱系統以提高熔化性能。建立了熔化傳熱過程的二維模型，分析了熔融前沿演化和動態溫度分布，研究了自然對流、翅片布局和Stefan 數的影響。通過響應面法以完全熔化時間為目標函數對翅片布置進行優化。結果表明，自然對流嚴重影響了相變儲熱裝置的熔化行為。Liu 等[18]為使相變儲熱系統的溫度更均勻、熔化速度更快，應用數值模擬方法優化了不均勻的樹狀肋片，研究填充角和中心角梯度在熔化過程中的作用，結果表明，不均勻的樹狀翅片顯著提高了儲能性能，因為熱傳導增強超過了自然對流抑制。與傳統的樹狀翅相比，通過儲熱系統中自然對流和熱傳導的協同改進，不均勻的樹狀翅片使溫度分布更均勻和熔化速度更快。Safari 等[19]將樹形翅片布置在水平放置相變儲熱（latent heat storage,LHS）單元的下部PCM 中，結果表明，與均勻的樹形翅片相比，非均勻設計將熔化過程縮短了9.25%。然而，非均勻翅片布局僅考慮下層PCM 的傳導增強，上層PCM 的熔化過程由于缺乏有效的傳熱路徑而減慢。Huang等[20]研究樹形翅片的LHS 單元均勻布置和梯度布置兩種方式。通過可視化實驗結合三維數值模擬，對LHS 單元的熔化∕凝固熱性能進行綜合分析。結果表明，樹狀翅片有利于熱從點到面的擴散，打破了傳統LHS 單元的傳熱滯后，從而加快了熔化∕凝固速度。在熔化過程中，梯度樹形翅片加強了熱能存儲（thermal energy storage, TES）單元下部的傳熱，延長了上部區域的對流傳熱持續時間，促進熔融后期自然對流和熱傳導的協同強化。與均勻翅片布局相比，梯度樹形翅片有效提高了熔化速度，將熔化持續時間縮短了9%。然而，梯度樹形翅片的非均勻傳輸路徑不利于以熱傳導為主的凝固傳熱，與均勻樹形相比，LHS 單元的溫度梯度增加，凝固時間延長了57.4%。

當前，研究主要從樹狀翅片結構、布局分析傳熱過程[21-22]，相變儲熱系統蓄釋熱整個傳熱過程深入的機理分析較少，部分研究者僅考慮將熱能儲存起來，忽略了釋熱過程的重要性，未深入分析自然對流與熱傳導在蓄釋熱過程中的協同作用。因此，本文基于二維梯度樹狀肋相變儲熱系統分析其傳熱特性，考慮到溫度分布的均勻性，采用場協同理論[23-25]分析溫度場和速度場的協同特性[26-28]，對相變材料的熔化和凝固過程的傳熱機理進行綜合分析，分析相變材料熔化溫度對蓄釋熱的影響，為相變材料的選擇提供理論依據。

1 模型的建立

1.1 物理模型

圖1為水平放置的三種相變儲熱裝置的二維模型。外管選擇鋁管，內直徑為124.0 mm，厚度為3.0 mm；內管選擇銅管，內直徑為45.0 mm，厚度為2.5 mm。內管與外管同軸，周向均勻布置六根相同的梯形翅片，肋片管材為銅，Case 1為均勻分布的六縱肋儲熱模型，Case 2 為在Case 1 模型基礎上增加了次肋的雪花型肋儲熱模型，Case 3 為在Case 2 基礎上采用梯度布置的樹狀肋儲熱模型。HTF 流經管內，PCM 填充在管外和殼體之間的空隙中。HTF 沿軸向方向溫度梯度變化與徑向方向的溫度梯度相比較小，傳熱特性主要與管子徑向方向有關，因此可以將三維傳熱模型簡化為二維的傳熱模型[16]。考慮到應用的可靠性及儲熱特性，PCM 選擇石蠟，HTF選擇熔融鹽，基礎材料的熱物性如表1所示。

圖1 三種相變儲熱模型及尺寸Fig.1 Three-phase change heat storage models and dimensions

表1 PCM及鋁和銅的熱物理性質［29］Table 1 Thermophysical properties of PCM，aluminum and copper［29］

1.2 數學模型

為簡化計算過程，對相變儲熱系統PCM 的熔化和凝固過程做如下假設：（1）傳熱問題是二維非穩態的；（2）液體PCM 被認為是不可壓縮的牛頓型流體，均質且各向同性；（3）儲熱系統初始溫度均勻且外壁是絕熱的；（4）液體PCM 的自然對流為層流[30]；（5）PCM 熔化過程中由于溫度變化引起密度變化使得浮升力作用下的傳熱服從Boussinesq 近似建模，在除浮力項外的所有項中密度均為常數；（6）PCM包括液相、固相和糊狀相；（7）液固界面處的兩相保持熱平衡；（8）HTF 以恒定質量流量和恒定的溫度進入管內。

在邊界條件方面，設置管內壁為等溫邊界條件[31]。管外及肋片與PCM 為流固耦合界面，儲熱罐外壁絕熱。因固液相變材料之間存在密度差，需要考慮重力引起的自然對流換熱，在Y軸方向開啟重力項，模型區域初始化定義，蓄熱過程相變材料及肋片管初始溫度為300.00 K，管內壁溫度設置為380.00 K[17]；釋熱過程相變材料及肋片管初始溫度為380.00 K，管內壁溫度設置為300.00 K。

1.3 數值分析方法

本模擬采用Fluent19.2 軟件，相變過程引入焓-孔隙率方法進行數值求解建模。采用二維焓-孔隙率法數值求解該問題，已有文獻指出該方法的結果與復雜的三維流體體積法的結果準確匹配[32]。采用SIMPLE 算法基于壓力-速度耦合的雙精度求解器，選擇壓力交錯選項PRESTO 用于自然對流，動量方程和能量方程的對流項采用二階迎風格式離散。為了更好地收斂，密度、壓力、動量和能量的欠松弛因子分別設置為1、0.3、0.5 和0.8。此外，當對應的連續性方程、動量方程及能量方程的最大殘差依次小于10-3、10-6、10-8，認為迭代是收斂的。

1.4 場協同理論分析傳熱過程

從對流換熱的場協同理論分析，對流傳熱系數是由整個對流區域中各種場參數相互作用的結果。對流傳熱的強度不僅取決于流體與固壁的溫差、流動速度、流體的熱物性和輸運性質，而且還取決于流體速度場與熱流場的協同程度。因此，引入場協同數

式中，Fc反映速度場和熱流場的協同程度（Fc≤1），當Fc=1時對流換熱的速度場和熱流場達到完全協同，此時，速度場與熱流場配合最好，即協同程度最高，同時也是換熱強度可能達到的最大值；-U為無量綱速度；?T為無量綱溫度梯度；β為流體流動方向與熱流傳遞方向的夾角，兩矢量的夾角盡可能小（β＜90°）或盡可能大（β＞90°）都能強化傳熱；速度矢量和熱流矢量的夾角的余弦總是趨于零（cosβ→0），使得Fc?1，因此通過提高速度場與溫度梯度場的協同可以強化換熱。

場協同數還可以表示成關于準則數的函數形式，其中Nusselt 數Nu為壁面上流體的無量綱溫度梯度，表征對流換熱的強弱；Grashof 數Gr判斷PCM熔化過程的流態；Prandtl數Pr是動量擴散與熱量擴散之比，表征流體的熱物理性能，具體表達式為式（2）～式（4）。

式中，h為對流傳熱系數，W∕(m2·K)；tw為翅片內管壁溫，380.00 K；t0為初始PCM 溫度，300.00 K；δ為PCM 填充沿徑向最大厚度，0.037 m，按照小空間自然對流換熱計算；運動黏度v=μ∕ρ，1.333×10-6m2∕s；a為熱擴散率，m2∕s；g為重力加速度，9.8 m2∕s；α為體積膨脹系數。本研究取翅片管壁溫與初始PCM 最大溫差為80℃，經計算Gr為1.8987×107，對于水平放置的儲熱系統，Gr在1.43×104～5.76×108為層流流動。

2 數值分析方法模型驗證

2.1 網格大小和時間步長獨立性驗證

利用結構化網格對二維截面進行網格劃分，如圖2所示，由四邊形網格和三角形網格組成。Case 2模型的熔化過程中，不同網格和時間步長PCM 液相率及溫度變化如圖3 所示。網格A、網格B、網格C節點數依次為12614、20480、34690，網格B 與網格C的溫度曲線幾乎重合，網格A 與網格B 的最大偏差為3.11%，網格B 與網格C 最大偏差僅為1.14%。時間步長0.05 s 和0.1 s 時的液相率具有很好的一致性，最大偏差僅為0.19%，與0.1 s 相比時間步長0.2 s 液相率最大偏差為1.84%。本研究選擇節點數20480，時間步長0.1 s。

圖2 計算域的網格圖Fig.2 Grid system of the computational domain

圖3 網格大小和時間步長獨立性驗證Fig.3 Verification of grid size and time step independence

2.2 模型準確性驗證

為驗證數值方法的有效性，將數值結果與Al-Abidi 等[33]的實驗結果進行比較，對三聯管換熱器（TTHX）的熔化進行了模型驗證。采用上述數值分析方法，選擇與文獻[33]相同的幾何模型和邊界條件。PCM 初始溫度為300.00 K。蓄熱開始，TTHX 內管壁由溫度為363.00 K 的HTF 加熱。如圖4 所示，內管、中管和外管的半徑分別為25.4 mm、75.0 mm、100.0 mm，當前模型PCM 平均溫度與實驗平均溫度比較，模擬數值與文獻的實驗數值吻合較好。實驗與模擬之間的最大相對誤差為2.11%，這是由于實驗測量不可避免的誤差，模擬過程進行了多條假設忽略了物性隨溫度的變化及固體下沉等因素，此誤差在允許的范圍內。因此，本文采用的數值模擬方法是有效的。

圖4 數值方法的驗證Fig.4 Validation of numerical method

3 結果分析

3.1 相變材料蓄釋熱過程云圖分析

3.1.1 蓄熱過程云圖分析 如圖5 所示，三個模型石蠟具有相同的熔化規律。初始階段，靠近管外壁及肋片附近區域的石蠟溫度先升高優先熔化，熱量的傳遞逐漸向外殼擴散，距離較近的區域仍處于低能區，這是由于石蠟的低熱導率使固態石蠟中的傳熱不利導致的，此時，熱傳導是傳熱的主要方式，對

圖5 三種模型石蠟熔化過程溫度及液相率云圖Fig.5 The cloud atlas of temperature and liquid fraction of paraffin melting process of three models

流換熱是微弱的。由于加熱管壁和相鄰石蠟之間存在較高的熱流量，焓的顯熱和潛熱被迅速吸收，因此，石蠟經歷了從固態—固液混合—液態的相變。隨著圍繞加熱管壁的液體石蠟的增加，在重力的作用下浮力驅動的自然對流使高溫低密度的液體石蠟產生向上的流動。高溫液體石蠟在殼體頂部區域的積累增加，與下部區域相比，吸熱量更高，此時自然對流換熱為傳熱的主要方式。隨著熔化的結束，液體石蠟在頂部區域的擁塞和較弱的流動性導致熱量傳遞減弱，熱傳導是換熱的主要方式。從云圖可知，Case 1 的石蠟熔化過程上部明顯比下部快，330 s 時仍有少量石蠟未熔化；Case 2、Case 3的熔化規律基本相同，這是由于在Case 1 的主肋上增設次肋進一步強化了傳熱，使熔化更均勻，290 s時石蠟全部熔化。可見，石蠟熔化過程是自然對流和熱傳導協同作用的結果。

3.1.2 釋熱過程云圖分析 如圖6 所示，在凝固初期，靠近肋片及管外壁的石蠟首先凝固，隨著非穩態過程的進行，釋熱的傳遞向外部擴散，石蠟放出熱量溫度降低，到達固相溫度開始凝固，100 s 石蠟未開始凝固，由初溫380.00 K 迅速下降，1200 s 時出現明顯的凝固。與熔化過程相比，石蠟在整個凝固過程所需時間較長，這是由于石蠟的熱導率低，凝固過程中熱傳導為主要的傳熱方式，對流換熱是極其微弱的[34]，因此，凝固過程溫度分布均勻。3600 s時Case 2、Case 3 石蠟全部凝固，Case 1 仍有部分石蠟未熔化。可見，在石蠟的凝固過程中熱傳導是換熱的主要方式。

圖6 三種模型石蠟凝固過程溫度及液相率云圖Fig.6 The cloud atlas of temperature and liquid fraction of paraffin solidification process of three models

3.2 相變材料蓄釋熱全過程傳熱特性分析

三種模型石蠟蓄釋熱全過程綜合分析如圖7～10 所示。由圖7 可知，三種模型的液相率在短時間內急劇升高，Case 1、Case 2、Case 3 的完全熔化時間依次為470、250、225 s，Case 3 所需時間最短，Case 1所需時間最長。三種模型釋熱過程所需時間較長，釋熱時間30 min，液相率依次為52.49%、10.10%、10.30%；釋熱時間50 min，液相率依次為37.43%、0%、0.71%。可見，Case 2 最先達到完全凝固，完全凝固時間為2918 s；其次是Case 3，完全凝固時間為3443 s；Case 1 所需時間最長。如圖8 所示，將石蠟的傳熱過程分為三個階段：蓄熱段、間歇段、釋熱段。在蓄熱段Case 3中石蠟平均溫度始終高于其他模型，在釋熱段石蠟的平均溫度均出現短時間急劇下降，釋熱時間100 s，三個模型石蠟的平均溫度由380.00 K 依次降為325.73、319.31、317.52 K；當石蠟平均溫度達到凝固溫度315.15 K，三個模型所需時間依次為254、169、145 s；達到凝固溫度后，石蠟的平均溫度與時間呈線性變化，Case 2、Case 3 幾乎重合，與Case 1 相比溫度變化更快。綜合分析蓄釋熱全過程，通過添加次肋強化了傳熱，肋片梯度分布對熔化過程有利，對凝固過程沒起到有利的作用。石蠟蓄釋熱過程具有熔化時間短、釋熱時間長的特點。一方面由于石蠟熔化過程伴隨著熱傳導與對流換熱的協同作用，釋熱過程對流換熱微弱以熱傳導為主；另一方面由于石蠟的熔化溫度315.15 K 更接近于300.00 K，距釋熱初溫380.00 K較遠。

圖7 三種模型蓄釋熱過程相變材料液相率變化Fig.7 Changes in the liquid fraction of phase change materials during heat storage and release process of three models

圖8 三種模型蓄釋熱過程相變材料平均溫度變化Fig.8 Changes in the average temperature of phase change materials during heat storage and release process of three models

三個模型Case 1、Case 2、Case 3 蓄釋熱過程石蠟平均速度瞬態變化如圖9 所示，蓄熱過程的速度變化呈先增大后減小的趨勢，初期和后期速度隨時間均呈線性變化，中期速度保持在一定范圍內波動并出現最大速度依次為：4.743、2.693、6.281 mm∕s，可見肋片的梯度分布有利于液體石蠟的擾動。釋熱過程初期，液體石蠟突然受到管內溫度為300.00 K 的流體冷卻，在開始的5 s內具有釋熱的最大速度依次為：11.069、4.389、8.254 mm∕s，這是由于液體石蠟突然受到邊界上的熱擾動產生的流動，而后速度下降直到釋熱約7 min 降低為零。綜合分析，石蠟的流動引起的自然對流換熱對蓄熱過程中期及釋熱過程初期起到重要的作用。如圖10所示，蓄釋熱過程熱通量均隨時間減小最終趨于零，這是由于非穩態傳熱石蠟與壁面的傳熱溫差逐漸減小，導致傳熱速度降低，當溫差為零時，蓄釋熱結束，無熱量傳遞。對比分析，Case3 的蓄釋熱熱通量最大。因此，梯度樹狀肋儲熱系統傳熱性能最優。

圖9 三種模型蓄釋熱過程相變材料速度變化Fig.9 Changes in the velocity of phase change materials during heat storage and release process of three models

圖10 三種模型蓄釋熱過程熱通量變化Fig.10 Changes in the heat flux density during heat storage and release process of three models

3.3 場協同理論分析蓄釋熱傳熱特性

如圖11、圖12所示，石蠟熔化過程中，Case 1 的場協同數基本穩定在0.02，且不隨時間發生較大的變化；Case 2 和Case 3 的場協同數隨蓄熱時間呈線性下降的趨勢，二者的場協同數大小相差不大；在時間τ＜150 s時，Case 2和Case 3的場協同數均大于Case 1，在50 s 時，Case 1、Case 2 和Case 3 的場協同數依次為：0.01964、0.05729、0.04869，而后Case 1 的場協同數略大于Case 2 和Case 3。石蠟的凝固過程中，Case 1 的場協同數變化不大，穩定在0.1 左右，Case 2 和Case 3 的場協同數隨釋熱的進行略有升高且均高于Case 1。這是由于與六縱肋片相比采用雪花型肋片和梯度樹狀肋片使空間溫度分布更均勻，影響到速度場的空間分布，從而加強了動量傳遞過程。說明在主肋上通過增設次肋可以提高流體速度場和溫度場的協同程度。

圖11 三種模型熔化過程石蠟場協同數Fc的瞬態變化Fig.11 Transient changes in synergy number Fc of paraffin field in melting process of three models

圖12 三種模型釋熱過程石蠟場協同數Fc的瞬態變化Fig.12 Transient changes in synergy number Fc of paraffin field during heat release of three models

3.4 相變材料熔化溫度對儲釋熱傳熱特性的影響

考慮到相變儲熱系統相變材料的蓄釋熱時間，以梯度樹狀肋儲熱系統（Case 3）為模型分析石蠟熔化溫度對蓄釋熱時間的影響。如圖13所示，以石蠟RT35、RT44、RT82 熔 點 溫 度 依 次 分 別 為315.00、340.00、360.00 K 進行分析[35-36]，相應的完全熔化時間為224、374、703 s；完全凝固時間為3439、1089、842 s。可見隨著石蠟熔化溫度的升高，完全熔化時間增長，完全凝固時間縮短。由圖14、圖15可知，蓄熱段低熔化溫度的石蠟溫升最快；釋熱過程開始，石蠟由初溫380.00 K 快速降低到達熔化溫度，而后需要較長的時間進行釋熱。釋熱過程初始階段石蠟熱通量較大，三個熔化溫度的最大熱通量依次為1662、1051、876 W∕m2，而后熱通量降低，最后趨于零釋熱結束，釋熱量隨熔化溫度的升高而減小。可見，熱通量越大溫度梯度越大，溫度曲線斜率越陡峭，熔化溫度為315.00、340.00 K 的石蠟先釋熱結束，熱通量的變化直接影響釋熱過程傳遞的熱量，在熱通量較大時釋熱速度較大，焓曲線先增加而后趨于平緩，且隨熔點溫度的升高而減小。因此，在選擇相變材料時要綜合考慮熔化溫度、蓄釋熱初溫和終溫及儲熱量的要求。

圖13 不同熔化溫度蓄釋熱過程石蠟液相率變化Fig.13 Changes in the liquid fraction of paraffin during heat storage and release at different melting temperatures

圖14 不同熔化溫度蓄釋熱過程石蠟平均溫度變化Fig.14 Changes in the average temperature of paraffin during heat storage and release at different melting temperatures

圖15 不同熔化溫度釋熱過程石蠟釋熱量變化Fig.15 Changes in the heat release capacity of paraffin at different melting temperatures

4 結 論

本文利用Fluent軟件對相變儲熱系統傳熱特性進行模擬，針對相變材料的熔化和凝固過程的傳熱機理進行綜合分析，對相變材料熔化溫度對蓄放熱時間的影響進行了研究，結合場協同理論對相變儲熱系統的溫度場和速度場的協同特性進行分析，得出以下結論。

（1）針對三種相變儲熱模型的蓄釋熱過程進行綜合分析。通過相變材料的液相率、平均溫度、平均速度及熱通量結合云圖分析蓄釋熱過程傳熱機理，石蠟熔化過程伴隨著熱傳導與自然對流的協同作用，凝固過程對流換熱微弱以熱傳導為主。綜合分析梯度樹狀肋儲熱系統傳熱性能最優。

（2）從場協同的角度分析了三種模型蓄釋熱傳熱特性，與六縱肋片相比，采用雪花型肋和梯度樹狀肋使空間溫度分布更均勻，影響到速度場的空間分布，從而加強了動量傳遞過程。說明在主肋上通過增設次肋可以提高流體速度場和溫度場的協同程度。

（3）相變材料的熔化溫度影響著蓄釋熱時間。石蠟熔化溫度分別為315.00、340.00、360.00 K，完全熔化時間依次為224、374、703 s；完全凝固時間依次為3439、1089、842 s。可見，隨著石蠟熔化溫度的升高，完全熔化時間增長，完全凝固時間縮短。釋熱過程開始，石蠟由初溫快速降到熔化溫度，而后需要較長的時間進行釋熱。熔化溫度為360.00 K的石蠟釋熱最慢，因此，在選擇相變材料時要考慮熔化溫度、蓄釋熱初溫和終溫及儲熱量的要求。