本文選自中國(guó)工程院院刊《Engineering》2021年第4期
作者：劉紅兵, 趙長(zhǎng)穎
來(lái)源：Effect of Radial Porosity Oscillation on the Thermal Performance of Packed Bed Latent Heat Storage[J].Engineering,2021,7(4):515-525.
編者按
太陽(yáng)能具有清潔、安全和取之不盡、用之不竭的巨大優(yōu)勢(shì)，已得到越來(lái)越廣泛的應(yīng)用，然而，太陽(yáng)能具有間歇性和波動(dòng)性，需要運(yùn)用儲(chǔ)能技術(shù)來(lái)解決其固有缺陷并提高能量存儲(chǔ)和利用效率。由于具有較高的儲(chǔ)熱能力和傳熱速率，堆積床相變儲(chǔ)熱被認(rèn)為是一種很有潛力的儲(chǔ)熱方法。
中國(guó)工程院院刊《Engineering》2021年第4期發(fā)表《堆積床相變儲(chǔ)熱系統(tǒng)中徑向孔隙率振蕩分布對(duì)熱性能的影響》。文章建立了一個(gè)基于球體實(shí)際堆積過(guò)程的三維堆積床相變儲(chǔ)熱模型，以描述徑向孔隙率的振蕩分布，并分析了其內(nèi)部的流動(dòng)和傳熱情況，揭示了徑向孔隙率振蕩分布的原因，探究了徑向孔隙率的分布對(duì)換熱流體的徑向相對(duì)速度分布以及相變材料的溫度分布和液相分?jǐn)?shù)分布的影響，并進(jìn)一步討論了無(wú)量綱參數(shù)直徑比、雷諾數(shù)、史蒂芬數(shù)對(duì)堆積床相變儲(chǔ)熱單元內(nèi)徑向流動(dòng)傳熱特性的影響。文章結(jié)果表明，不同的直徑比對(duì)應(yīng)于不同的徑向孔隙率分布；隨著直徑比的增加，在壁面附近區(qū)域換熱流體速度顯著變化，而中心區(qū)域換熱流體速度的不均勻性將減??；雷諾數(shù)和史蒂芬數(shù)對(duì)換熱流體的相對(duì)速度分布有輕微影響，而更高的雷諾數(shù)可導(dǎo)致速度成比例地提高，史蒂芬數(shù)的增加可加快堆積床相變儲(chǔ)熱系統(tǒng)的儲(chǔ)熱過(guò)程。
一、引言
太陽(yáng)能由于具有清潔、安全和取之不盡、用之不竭的巨大優(yōu)勢(shì)而得到了越來(lái)越廣泛的應(yīng)用。根據(jù) REN21的報(bào)告，聚光太陽(yáng)能發(fā)電（CSP）設(shè)施和太陽(yáng)能光伏發(fā)電設(shè)施（PV）的總裝機(jī)容量從2018年的512 GW 增長(zhǎng)至2019年的627 GW，占可再生能源增長(zhǎng)總量的 60.21%。
然而，太陽(yáng)能具有間歇性和波動(dòng)性，需要儲(chǔ)能技術(shù)來(lái)解決其固有缺陷并提高能量存儲(chǔ)和利用效率。在太陽(yáng)能的熱利用方面，相變儲(chǔ)熱擁有更高的儲(chǔ)能密度并在相變過(guò)程中保持幾乎恒定的溫度，該特性引起了廣泛的關(guān)注。相變儲(chǔ)熱的主要限制因素是相變材料較低的導(dǎo)熱系數(shù)（0.2~0.8 W·m?1·K?1），因此需要采取有效措施提升其傳熱速率。與殼管式單元相比，對(duì)于相變蓄熱單元，堆積床是一種更加高效的熱能存儲(chǔ)方式。據(jù)Li等研究，堆積床相變儲(chǔ)熱系統(tǒng)的儲(chǔ)熱速率是管殼式儲(chǔ)熱系統(tǒng)的1.8~3.2倍。
在堆積床相變儲(chǔ)熱系統(tǒng)中，由于換熱流體直接流過(guò)相變膠囊，因此相變膠囊的分布情況會(huì)極大地影響流場(chǎng)的分布。對(duì)于球形堆積，最密堆積方式為面心立方堆積和六方最密堆積，其堆積密度達(dá)到74.05%。但是，當(dāng)球體隨機(jī)分布在容器中時(shí)，由于壁面的緣故，球體無(wú)法達(dá)到最密排列結(jié)構(gòu)，這種現(xiàn)象稱為壁面效應(yīng)。Mueller的研究表明，堆積床內(nèi)部孔隙率沿著徑向方向呈振蕩分布。此外，徑向孔隙的振蕩分布會(huì)導(dǎo)致流場(chǎng)不均勻，這將進(jìn)一步影響堆積床內(nèi)部的傳熱效果。由于容器為圓柱形，壁面對(duì)膠囊的平均填充密度的影響隨著容器和球體的直徑比（D/dp）的減小而顯著增加。當(dāng)使用堆積床進(jìn)行相變儲(chǔ)熱時(shí)，如果想要達(dá)到較高的儲(chǔ)熱能力則要求系統(tǒng)具有較小的儲(chǔ)罐尺寸和直徑比。因此，在堆積床相變儲(chǔ)熱系統(tǒng)的數(shù)值研究中應(yīng)考慮徑向孔隙率的變化，特別是對(duì)于小直徑比（D/dp<10）的情況。
針對(duì)堆積床相變儲(chǔ)熱的計(jì)算模型，根據(jù)是否需要描述具體的堆積結(jié)構(gòu)，可以將模型分為兩類(lèi)：基于表征體元尺度（REV）的模型和基于孔隙尺度的模型。在基于表征體元尺度的模型中，主要通過(guò)控制方程計(jì)算換熱流體與相變材料的熱物性變化。Benmansour等采用雙方程模型來(lái)表示換熱流體與相變材料的能量方程，而方程中采用的孔隙率是一個(gè)常數(shù)。因此，該研究沒(méi)有考慮徑向孔隙率的變化。隨后的研究中，他們用一個(gè)指數(shù)表達(dá)式替代了能量方程中屬于常數(shù)項(xiàng)的孔隙率。在該表達(dá)式中，徑向孔隙率φ(r)是與壁面相隔距離的函數(shù)。此外，為了計(jì)算相變膠囊內(nèi)部的溫度梯度，有研究提出了Dispersion-concentric (DC)模型。Karthikeyan和Velraj比較了雙方程模型和DC模型；其結(jié)果表明，考慮到相變膠囊內(nèi)部的內(nèi)部傳熱，DC模型的結(jié)果與實(shí)驗(yàn)結(jié)果更加吻合。但是，其所采用的數(shù)值方法是REV方法，而小直徑比（D/dp<10）的堆積床不能視為連續(xù)的多孔介質(zhì)。由于網(wǎng)格無(wú)法精確描述孔隙和實(shí)體骨架，因此無(wú)法滿足REV尺度模型的要求。
在孔隙尺度模型中，基于對(duì)相變膠囊的具體分布的描述，可以計(jì)算出換熱流體的流量以及換熱流體和相變膠囊之間的傳熱量。Xia等開(kāi)發(fā)了一種高效的堆積床模型，將三維堆積床轉(zhuǎn)化為二維堆積床。相應(yīng)的轉(zhuǎn)化標(biāo)準(zhǔn)是確保兩個(gè)堆積床具有相同的孔隙率。但是，在轉(zhuǎn)換過(guò)程中該模型未考慮徑向孔隙率的分布。此外，該模型中的相變膠囊并不相互接觸，并且忽略了膠囊之間的熱傳導(dǎo)。由于轉(zhuǎn)化過(guò)程中相變膠囊之間的流道發(fā)生了變化，因此其流場(chǎng)與真實(shí)流場(chǎng)有差異。
綜上所述，本文根據(jù)球體的實(shí)際堆積過(guò)程，構(gòu)建了考慮徑向孔隙率的振蕩的三維堆積床相變儲(chǔ)熱模型，揭示了徑向孔隙率振蕩分布的原因，探究了徑向孔隙率的分布對(duì)換熱流體的徑向相對(duì)速度分布以及相變材料的溫度分布和液相分?jǐn)?shù)分布的影響，并進(jìn)一步討論了無(wú)量綱參數(shù)直徑比、雷諾數(shù)、史蒂芬數(shù)對(duì)堆積床相變儲(chǔ)熱單元內(nèi)徑向流動(dòng)傳熱特性的影響。
二、三維堆積床相變儲(chǔ)熱模型的建立
（一）物理模型的建立
在研究三維堆積床時(shí)，面臨的主要困難是建立有效的球體填充模型。球體堆積過(guò)程中，小球之間存在碰撞和摩擦作用，因此需要計(jì)算小球在重力、摩擦力以及相互作用力的作用下的運(yùn)動(dòng)軌跡，最終確定穩(wěn)定狀態(tài)下所有小球的位置。
建立堆積床的現(xiàn)有方法包括DEM方法、Blender法和蒙特卡洛法。在本研究中，采用開(kāi)源軟件Blender建立了三維堆積床模型。 Blender軟件搭載了Bullet Physics Library，可以有效地計(jì)算小球在重力作用下的物理過(guò)程。在建模過(guò)程中，首先建立了一個(gè)圓柱形容器。然后，不斷在容器上方隨機(jī)產(chǎn)生小球，使其掉落，如圖1（a）所示。在堆積過(guò)程中計(jì)算每個(gè)球的位置。如果所有球的位置都達(dá)到穩(wěn)定狀態(tài)，則可以認(rèn)為堆積過(guò)程達(dá)到了穩(wěn)定狀態(tài)。
孔隙率是空隙體積與總體積的比率，而表面孔隙率是空隙面積與表面總面積的比率，如等式（1）~（3）所示。在這項(xiàng)研究中討論的徑向孔隙率是表面孔隙率。為了計(jì)算徑向孔隙率，沿徑向方向截取了一系列圓柱表面。對(duì)于每個(gè)圓柱表面，采用射線投射算法來(lái)確定球體內(nèi)部區(qū)域。然后將所有球體內(nèi)部面積相加，并用球體的總面積除以圓柱表面的面積，如圖1（b）所示。
圖1 （a）球的堆積過(guò)程；（b）沿徑向生成次表面的示意圖
在圖2中，在Blender中建立了直徑比為7.99的模型，將徑向孔隙率與Mueller的結(jié)果進(jìn)行了比較。由圖2 可知，兩條曲線均呈振蕩分布且吻合較好，從而驗(yàn)證了模型的可靠性。此外，徑向孔隙率的振蕩分布表明堆積床不能處理成均勻多孔介質(zhì)，這也進(jìn)一步證明了建立三維模型的必要性。
圖2 直徑比為7.99時(shí)徑向孔隙率分布的比較
在堆積床中，球之間的接觸類(lèi)型為點(diǎn)接觸，這不利于網(wǎng)格的生成。當(dāng)前有四種方法處理接觸點(diǎn)問(wèn)題：重疊法、間隙法、壓平法和橋梁法。在具有熱傳導(dǎo)的情況下，橋梁法可以達(dá)到最佳效果。因此，直徑為dp/10的小圓柱體作為熱橋被插入到每?jī)蓚€(gè)接觸的球體之間，這會(huì)導(dǎo)致堆積床的平均孔隙度變化0.05%~0.08%。此外，球體與壁面接觸點(diǎn)的處理采用了間隙法，在球體和容器壁之間形成一個(gè)1 mm的間隙。
如圖3所示，本文建立了三維堆積床相變儲(chǔ)熱模型。在該模型中，將205個(gè)相變膠囊裝入一個(gè)內(nèi)徑為240 mm 的圓柱形罐中，其中相變膠囊的直徑為48 mm；計(jì)算得到堆積床的直徑比為5。選擇三元碳酸鹽Li2CO3- K2CO3-Na2CO3（32 wt%-35 wt%-33 wt%）作為相變材料，其熱物性如表1所示。堆積床被放置在中間部分，同時(shí)在入口和出口的兩側(cè)添加延伸部分以消除進(jìn)出口段的影響。
表1 三元碳酸鹽相變材料熱物性
Ts: solid state temperature; Tl : liquid state temperature; ΔH : specific enthalpy; λ: thermal conductivity; cp: Specific heat capacity; ρ : density; s: solid state; l: liquid state.
圖3 三維堆積床相變儲(chǔ)熱系統(tǒng)模型示意圖
（二）數(shù)學(xué)模型
1. 模型假設(shè)
（1）忽略重力對(duì)換熱流體的影響；
（2）忽略相變膠囊中相變材料的對(duì)流作用；
（3）忽略換熱流體與相變膠囊的輻射傳熱。
2. 換熱流體的控制方程
在堆積床相變儲(chǔ)熱系統(tǒng)中，基于相變膠囊直徑的雷諾數(shù)定義為：
式中，ρ是密度；μ是黏度；uin是入口速度。
由于堆積床內(nèi)流道復(fù)雜，內(nèi)部流動(dòng)通常為湍流。考慮到相變膠囊的換熱表面為球面，并且為了增強(qiáng)對(duì)壁面的處理，換熱流體的流動(dòng)使用RNG k-ε 湍流模型進(jìn)行計(jì)算。因此堆積床相變儲(chǔ)熱系統(tǒng)的控制方程在直角坐標(biāo)系下可以寫(xiě)成如下形式。
連續(xù)性方程：
動(dòng)量方程：
式中，xi 和xj 分別是在i 和j 方向上的距離；ui 和uj 分別是換熱流體在i和j方向上的速度；ui′ 和uj′分別是換熱流體在i 和j 方向上的波動(dòng)速度；E是能量；(τij)eff是有效的雷諾張量；k是湍流動(dòng)能；ε是湍流耗散率；σT 是能量方程式中的普朗特?cái)?shù)；σk 是湍流動(dòng)能方程中的普朗特?cái)?shù)；σε是湍流耗散率方程式中的普朗特?cái)?shù)；C1ε和C2ε是兩個(gè)模型常數(shù)。湍流黏度μt 、湍流動(dòng)能的產(chǎn)生項(xiàng)Gk以及附加項(xiàng)Rε由下式給出：
式中，Cμ是模型常數(shù)；γ是RNG k-ε湍流模型系數(shù)；γ0是初始狀態(tài)下的RNG k-ε湍流模型系數(shù)；S是流體的平均應(yīng)變率；Sij 是變形張量。
上述方程式中的常數(shù)為：σT = 1.0，σk =σε = 0.72，Cμ = 0.0845，C1ε = 1.42，C2ε =1.68，γ0 = 4.38。
3. 相變材料的控制方程
相變膠囊中相變材料的自然對(duì)流作用較弱，故而忽略，因此能量方程可以表示為：
式中，Hpcm代表相變材料的焓值，定義為顯熱和潛熱之和；? 是梯度算子。
式中，H0,pcm代表相變材料在初始溫度下的焓值；ΔH代表相變材料儲(chǔ)熱過(guò)程中釋放的相變潛熱，定義為ΔH = βL，L代表相變材料的潛熱值，β代表相變材料的液相分?jǐn)?shù)，β在不同的溫區(qū)下表示為：
式中，Tpcm, T0, Ts, Tl分別代表相變材料的溫度、初始溫度、開(kāi)始熔化溫度和終止熔化溫度。
另外，儲(chǔ)熱過(guò)程中相變材料的史蒂芬數(shù)定義為顯熱量與潛熱量的比值：
式中，Tin是指入口換熱流體的溫度。
4. 邊界條件與初始條件
儲(chǔ)熱過(guò)程開(kāi)始前，堆積床相變儲(chǔ)熱系統(tǒng)保持在初始溫度且有Tf =Tpcm =T0。儲(chǔ)熱開(kāi)始時(shí)，入口處的換熱流體以恒定的溫度Tin和質(zhì)量流量qin進(jìn)入系統(tǒng)中，設(shè)出口處換熱流體在垂直于出口方向的熱流為零。因此，堆積床相變儲(chǔ)熱系統(tǒng)中邊界條件如下。
進(jìn)口處：
出口處：
式中，z 是高度。
同時(shí)為了后續(xù)更好地分析徑向孔隙率對(duì)換熱流體徑向特性的影響，儲(chǔ)熱罐的壁面設(shè)置為絕熱邊界條件：
（三）堆積床相變儲(chǔ)熱系統(tǒng)的儲(chǔ)熱
堆積床相變儲(chǔ)熱系統(tǒng)的總儲(chǔ)熱量包括相變材料、相變膠囊殼體以及儲(chǔ)熱罐壁面存儲(chǔ)的熱量。其中，相變材料儲(chǔ)熱量包括固態(tài)顯熱、液態(tài)顯熱和相變潛熱，而相變膠囊殼體以及儲(chǔ)熱罐壁面存儲(chǔ)的熱量均為顯熱。系統(tǒng)的儲(chǔ)熱時(shí)間定義為靠近出口處相變膠囊溫度達(dá)到低于入口溫度1 K所用的時(shí)間。因此，堆積床相變儲(chǔ)熱系統(tǒng)的總儲(chǔ)熱量Qtotal可以表示為：
式中，Qpcm、Qshell和Qtank分別是相變材料、相變材料膠囊殼和儲(chǔ)罐中的熱量存儲(chǔ)；mpcm、mshell和mtank分別是相變材料的質(zhì)量、相變材料膠囊的外殼和儲(chǔ)罐的質(zhì)量；Cp,s和Cp,l是固態(tài)和液態(tài)相變材料的比熱容；Tm、T0和Tin分別是熔化溫度、初始溫度和入口溫度；Cp,shell和Cp,tank分別是外殼和儲(chǔ)罐的熱容量。
在儲(chǔ)熱過(guò)程中，堆積床相變儲(chǔ)熱系統(tǒng)的實(shí)際儲(chǔ)熱量 Qstored可以通過(guò)熱力學(xué)第一定律計(jì)算得到：
式中，mf是換熱流體的質(zhì)量流速；Cpf,in和Cpf,out分別是換熱流體在入口和出口的比熱容；Tf,in和Tf,out分別是入口和出口處換熱流體的溫度。
儲(chǔ)熱過(guò)程中的平均儲(chǔ)熱功率為：
式中，τcharge是儲(chǔ)熱時(shí)間。
（四）數(shù)值計(jì)算方法
由于有數(shù)百個(gè)相變膠囊隨機(jī)填充在堆積床中，因此采用非結(jié)構(gòu)化網(wǎng)格對(duì)該三維結(jié)構(gòu)進(jìn)行網(wǎng)格劃分。首先利用shrink-wrap方法對(duì)堆積床球體表面進(jìn)行包面處理，然后使用多面體網(wǎng)格進(jìn)行網(wǎng)格的填充。隨后，確定網(wǎng)格尺寸對(duì)填充時(shí)間和堆積床中壓降的影響，并選擇dp/18 網(wǎng)格尺寸進(jìn)行網(wǎng)格生成，如圖4所示。在計(jì)算過(guò)程中，利用SIMPLE算法實(shí)現(xiàn)壓力場(chǎng)與速度場(chǎng)的耦合迭代，離散方程采用二階迎風(fēng)格式，瞬態(tài)方程的求解則是使用了二階隱式格式。
圖4 網(wǎng)格無(wú)關(guān)性驗(yàn)證
Li等搭建了堆積床相變儲(chǔ)熱實(shí)驗(yàn)系統(tǒng)，該系統(tǒng)由385個(gè)相變膠囊隨機(jī)堆積而成。系統(tǒng)中儲(chǔ)熱罐內(nèi)徑為 260 mm，相變膠囊外徑為34 mm，所以儲(chǔ)熱罐與相變膠囊的直徑比為7.65。相變膠囊中的相變材料選取三元碳酸鹽 Li2CO3-K2CO3-Na2CO3（32 wt%-35 wt%-33 wt%），其熱物性如表1所示。該系統(tǒng)測(cè)試了三個(gè)不同高度方向上相變膠囊內(nèi)部相變材料在儲(chǔ)熱過(guò)程中的溫度變化，同時(shí)可以驗(yàn)證本文建立的三維堆積床相變儲(chǔ)熱系統(tǒng)模型。如圖5（a）所示，相變膠囊沿高度的位置為z/H = 0.25、 0.5、0.75。在裝料過(guò)程中，質(zhì)量流量qin為260 kg·h?1，而堆積床的初始溫度T0為325 ℃，換熱流體的入口溫度 Tin為465 ℃。
如圖5（b）所示，三個(gè)位置處的模擬數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)均吻合較好，其中z/H = 0.5處的結(jié)果最佳。儲(chǔ)熱過(guò)程剛開(kāi)始時(shí)，模擬數(shù)據(jù)快于實(shí)驗(yàn)數(shù)據(jù)；但是相變材料完成相變過(guò)程后，模擬數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)基本同步，并且模擬與實(shí)驗(yàn)幾乎同時(shí)完成儲(chǔ)熱過(guò)程。一方面，將換熱流體從初始溫度（如325 ℃）加熱到入口溫度（如465 ℃）的加熱過(guò)程較長(zhǎng)，而在模擬中這一過(guò)程是在很短時(shí)間內(nèi)完成的。另一方面，當(dāng)儲(chǔ)熱過(guò)程進(jìn)行到液態(tài)儲(chǔ)熱階段時(shí)，相變材料已經(jīng)完成了相變儲(chǔ)熱過(guò)程，初始狀態(tài)下入口溫度變化的不同而造成的儲(chǔ)熱量的差異相對(duì)于系統(tǒng)整體的儲(chǔ)熱量而言已經(jīng)很小，所以模擬數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)越往后吻合度越高。總體而言，模擬與實(shí)驗(yàn)結(jié)果的對(duì)比驗(yàn)證了三維堆積床相變儲(chǔ)熱模型的可靠性與有效性。
圖5（a）熱電偶布置示意圖；（b）模擬結(jié)果與實(shí)驗(yàn)數(shù)據(jù)對(duì)比
三、結(jié)果與討論
（一）流動(dòng)性能與熱性能分析
本節(jié)中，堆積床的初始溫度T0為325 ℃，直徑比為 5，HTF的入口溫度為465 ℃，質(zhì)量流量為260 kg· h?1 。因?yàn)槎逊e床的平均孔隙率為0.437，計(jì)算得到相應(yīng)的雷諾數(shù)Rep為4998.72，可認(rèn)為堆積床內(nèi)的流動(dòng)為湍流。
1. 徑向孔隙率振蕩分布的原因
對(duì)于直徑比為5的堆積床相變儲(chǔ)熱單元，其徑向孔隙率分布如圖6所示。為了更直觀地說(shuō)明，截取了堆積床內(nèi)部不同徑向位置的圓柱面，這些圓柱截面距離儲(chǔ)熱罐壁面的距離分別為0、0.25dp、0.5dp、0.75dp、dp、 1.25dp、1.5dp、1.75dp、2dp和2.25dp。圓柱截面上藍(lán)色區(qū)域代表?yè)Q熱流體，白色區(qū)域則為相變材料。根據(jù)面孔隙率的定義可知，距離壁面0、dp和2dp的圓柱截面孔隙率高于其他截面的孔隙率，而這與相變膠囊的堆積過(guò)程以及膠囊的形狀有關(guān)。在堆積過(guò)程中，相變膠囊往往傾向于先填充靠近墻壁的空間，然后再填充儲(chǔ)熱罐中心的區(qū)域。由于相變膠囊的外殼是球形，所以壁面處相變膠囊與儲(chǔ)熱罐壁面只有點(diǎn)接觸，該位置的孔隙率則接近 1，這也意味著圓柱截面上幾乎全是換熱流體區(qū)域。在距離壁面0.5dp處，圓柱截面剛好穿過(guò)相變膠囊的中心，這意味著相變材料區(qū)域?qū)⒄紦?jù)圓柱截面的大部分區(qū)域，從而導(dǎo)致孔隙率降低到了0.15。同樣地，在距離壁面dp 的圓柱截面上，靠近壁面的相變膠囊與內(nèi)部相變膠囊的接觸方式也是點(diǎn)接觸，所以在接觸位置附近能夠形成較大的間隙，導(dǎo)致孔隙率增大至0.67。因此，徑向圓柱截面上的孔隙率表現(xiàn)為高低交錯(cuò)的分布，即徑向孔隙率呈振蕩分布。
圖6 徑向孔隙率分布及徑向圓柱截面示意圖
2. 流速分布
類(lèi)似于堆積床中徑向孔隙率的定義，徑向流速為不同徑向位置圓柱截面上換熱流體的平均流速，定義為：
式中，Ωf代表圓柱截面上換熱流體的區(qū)域；Af代表圓柱截面上換熱流體區(qū)域的面積。
相對(duì)流速U定義為徑向流速u(mài)(r)與入口流速u(mài)in之比：
圖7展示了堆積床內(nèi)換熱流體相對(duì)速度沿徑向方向的分布，可以發(fā)現(xiàn)相對(duì)速度也呈現(xiàn)出振蕩的分布。盡管堆積床在壁面處的孔隙率接近1，但由于換熱流體與壁面之間存在黏性作用，壁面處的流速為0。沿著遠(yuǎn)離壁面的方向，換熱流體的流速急劇增大并呈現(xiàn)出與徑向孔隙率一致的分布。換熱流體進(jìn)入堆積床后從相變膠囊之間的空隙流過(guò)，空隙越大的徑向位置對(duì)換熱流體的阻力就越小。所以孔隙率大的位置有更多的換熱流體流過(guò)，相應(yīng)的流速也更大。為了便于對(duì)相對(duì)流速分布進(jìn)行分析與比較，我們將堆積床沿徑向方向分為兩個(gè)區(qū)域：近壁面區(qū)域和中心區(qū)域。其中近壁面區(qū)域?yàn)閺谋诿娴竭h(yuǎn)離壁面0.5dp的區(qū)域，流速在這里發(fā)生了劇烈的變化；而中心區(qū)域則是從遠(yuǎn)離壁面0.5dp到罐體中心的區(qū)域，流速在這里的變化則更為平緩。
圖7 徑向方向相對(duì)速度分布
近壁面區(qū)域流速變化劇烈，所以評(píng)判該區(qū)域最重要的參數(shù)便是壁面附近的最大相對(duì)速度。在本算例中，壁面附近最大的相對(duì)速度為3.49。根據(jù)連續(xù)性方程可知，在相同的流量下，流速的變化與流道面積變化成反比 [見(jiàn)式（31）]。因此當(dāng)流體在堆積床中流動(dòng)時(shí)，如果流道突然收縮，流速會(huì)相應(yīng)地增大。圖8展示了堆積床中不同位置截面上的流速分布。能夠明顯觀察到，壁面附近換熱流體的流速在經(jīng)歷流道收縮時(shí)大幅增大，這也就是壁面附近流速變化劇烈的原因。
圖8 堆積床內(nèi)流速分布
中心區(qū)域換熱流體的相對(duì)流速呈現(xiàn)出更加平緩的分布，和徑向孔隙率的分布一致。為了衡量換熱流體流速的不均勻性，我們提出使用相對(duì)流速的標(biāo)準(zhǔn)差進(jìn)行評(píng)估，如式（32）所示。標(biāo)準(zhǔn)差越大意味著徑向位置上流速的分布差異越大，即流速的分布越不均勻。對(duì)于該直徑比為5的堆積床，當(dāng)雷諾數(shù)Rep為4998.72時(shí)，相對(duì)流速的標(biāo)準(zhǔn)差為0.39。
式中，N為截面的個(gè)數(shù)；ui為第i個(gè)截面的相對(duì)速度；u為所有截面相對(duì)流速的平均值。
3. 徑向溫度和液相分?jǐn)?shù)分布
相變材料的徑向溫度和徑向液相分?jǐn)?shù)分別定義如下：
式中，Ωpcm代表圓柱截面上相變材料的區(qū)域；Apcm代表圓柱截面上相變材料區(qū)域的面積。
圖9（a）和（b）分別展示了不同時(shí)刻下相變材料徑向溫度分布和徑向液相分?jǐn)?shù)分布?？梢钥吹剑煌恢锰幭嘧儾牧系臏囟群鸵合喾?jǐn)?shù)均不一致，意味著不同位置處相變材料的儲(chǔ)熱進(jìn)程不一致。當(dāng)t = 200 s時(shí)，除壁面附近外，其他徑向位置上相變材料溫度還未達(dá)到熔點(diǎn)，因此液相分?jǐn)?shù)保持為0，此時(shí)相變材料還處于固態(tài)儲(chǔ)熱階段。當(dāng)t = 1200 s時(shí)，徑向位置上相變材料的溫度處于相變溫度范圍內(nèi)，因此液相分?jǐn)?shù)出現(xiàn)了很大的波動(dòng)。值得注意的是，由于相變材料的相變溫度范圍為 395.1~413 ℃，根據(jù)焓值法可知，每1 ℃的升溫均會(huì)引起液相分?jǐn)?shù)增加0.056。所以，雖然各位置處相變材料在相變過(guò)程時(shí)溫度差別不大，但液相分?jǐn)?shù)會(huì)發(fā)生較大幅度的變化。當(dāng)t = 2400 s時(shí)，各個(gè)徑向位置處相變材料的溫度基本都超過(guò)了相變溫度范圍上限，所以此時(shí)液相分?jǐn)?shù)幾乎都為1，大部分相變材料處于液態(tài)儲(chǔ)熱階段。
圖9 相變材料不同時(shí)刻的徑向分布。（a）溫度；（b）液相分?jǐn)?shù)
圖10刻畫(huà)了t = 600 s時(shí)距離壁面為0.5 dp、dp、1.5dp 和2dp的圓柱截面上相變材料的液相分?jǐn)?shù)分布。首先可以發(fā)現(xiàn)，靠近入口處的相變膠囊中的相變材料熔化更快，并且熔化從換熱流體的流動(dòng)側(cè)開(kāi)始。其次，距離壁面為dp和2dp的圓柱截面上相變材料熔化更快，與圖9 （b）中液相分?jǐn)?shù)的分布情況相符。一方面，由于0.5dp 和1.5dp對(duì)應(yīng)的圓柱截面孔隙率低，圓柱截面上換熱流體的流速相對(duì)較低，所以對(duì)流換熱作用更弱；另一方面， 0.5dp和1.5dp對(duì)應(yīng)的圓柱截面上相變材料占比大，熔化需要的熱量更多，所以儲(chǔ)熱進(jìn)程更為緩慢。
圖10 t = 600 s時(shí)徑向圓柱截面上相變材料的液相分?jǐn)?shù)
（二）直徑比對(duì)堆積床流動(dòng)換熱的影響
本節(jié)研究了具有不同直徑比（4, 5, 6）的堆積床相變儲(chǔ)熱系統(tǒng)的流動(dòng)特性和傳熱儲(chǔ)熱性能。入口溫度Tin為 465 ℃，質(zhì)量流率為260 kg· h?1 儲(chǔ)熱罐的內(nèi)徑為240 mm，高度為500 mm。為了便于分析比較，三個(gè)堆積床相變儲(chǔ)熱系統(tǒng)的儲(chǔ)熱量在設(shè)計(jì)時(shí)基本一致，并基于此確定了三個(gè)系統(tǒng)中所需的相變膠囊個(gè)數(shù)。三個(gè)不同直徑比堆積床內(nèi)相變膠囊的尺寸、個(gè)數(shù)和孔隙率等參數(shù)總結(jié)在表2 中。圖11展示了不同直徑比堆積床的徑向孔隙率分布，通過(guò)比較發(fā)現(xiàn)：具有不同直徑比的堆積床對(duì)應(yīng)于不同的徑向孔隙率分布。從壁面到距離壁面1.5dp的范圍內(nèi)，三個(gè)堆積床的徑向孔隙率的分布基本一致，但是越靠近罐體中心，徑向孔隙率的差別越大。對(duì)于直徑比為4的堆積床，容器的直徑正好能夠容納4個(gè)相變膠囊，所以相變膠囊在罐體中心是以點(diǎn)接觸或者不接觸的形式排布，導(dǎo)致罐體中心區(qū)域的孔隙率很大。對(duì)于直徑比為5的堆積床，容器的直徑正好能夠容納5個(gè)相變膠囊，所以罐體中心被相變膠囊所占據(jù)，使得罐體中心的孔隙率幾乎為0。對(duì)于直徑比為6的堆積床，容器的直徑正好能夠容納6個(gè)相變膠囊，理論上在罐體中心會(huì)與直徑比為4的堆積床一樣具有很高的孔隙率。但是直徑比為6時(shí)壁面效應(yīng)的影響減弱，罐體中心相變膠囊的排列具有一定的隨機(jī)性，如圖11所示。
表2 不同直徑比的堆積床相變儲(chǔ)熱系統(tǒng)的參數(shù)
圖11 不同直徑比堆積床的徑向孔隙率分布
圖12（a）展示了具有不同直徑比的堆積床的徑向相對(duì)速度的分布，與徑向孔隙率的分布基本一致。在近壁面區(qū)域，直徑比為4、5、6的堆積床最大相對(duì)速度依次增大，分別為2.87、3.49和3.99。相變膠囊的直徑隨直徑比的增大而減小，使得相變膠囊與壁面間的空隙變狹窄。由式（31）可知，當(dāng)熱空氣流經(jīng)狹窄的流道時(shí)，流速會(huì)增大。在中心區(qū)域，三個(gè)堆積床內(nèi)相對(duì)速度的標(biāo)準(zhǔn)差分別為0.53、0.39和0.27，這意味著隨著直徑比的增大，堆積床中心區(qū)域流速變得更加均勻。
圖12（b）展示了t = 1200 s時(shí)具有不同直徑比的堆積床的相變材料徑向液相分?jǐn)?shù)的分布，能夠觀察到：直徑比為4時(shí)堆積床近壁面處與罐體中心處相變材料具有最快的儲(chǔ)熱速度；而對(duì)于直徑比為5和6的堆積床，近壁面處的相變材料的儲(chǔ)熱速度明顯快于罐體中心區(qū)域的相變材料，這與徑向孔隙率的分布是一致的。
圖12 不同直徑比相變堆積床LHS系統(tǒng)的徑向分布。（a）相對(duì)流速；（b） t = 1200 s時(shí)不同直徑比堆積床內(nèi)相變材料徑向液相分?jǐn)?shù)分布
圖13展示了不同直徑比堆積床相變儲(chǔ)熱系統(tǒng)的儲(chǔ)熱量與平均儲(chǔ)熱功率隨時(shí)間的變化，直徑比分別為4、5 和6時(shí)系統(tǒng)的儲(chǔ)熱時(shí)間分別為121.03 min、89.73 min和 79.80 min。由于三個(gè)儲(chǔ)熱系統(tǒng)在設(shè)計(jì)時(shí)的儲(chǔ)熱量基本一致，所以儲(chǔ)熱系統(tǒng)的平均儲(chǔ)熱功率隨著直徑比的增大而增大。直徑比為4的堆積床中相變膠囊的直徑為60 mm，而直徑比為5和6的堆積床中相變膠囊的直徑分別為 48 mm和40 mm。小的直徑對(duì)應(yīng)更短的傳熱距離，因此儲(chǔ)熱進(jìn)程能夠加速。但是相變膠囊直徑越小，堆積床的結(jié)構(gòu)越緊密，這會(huì)引起堆積床內(nèi)壓降的急劇增大。計(jì)算結(jié)果顯示，直徑比分別為4、5和6的堆積床內(nèi)的壓降分別為142.97 Pa、213.34 Pa和257.13 Pa，這意味著維持系統(tǒng)運(yùn)行所需的泵功依次增大。因此，當(dāng)考慮將堆積床相變儲(chǔ)熱系統(tǒng)應(yīng)用于實(shí)際場(chǎng)景時(shí)，需要綜合考慮儲(chǔ)熱系統(tǒng)的儲(chǔ)熱時(shí)間、儲(chǔ)熱功率以及內(nèi)部壓降的大小等因素。
圖13 不同直徑比堆積床相變儲(chǔ)熱系統(tǒng)的儲(chǔ)熱量與平均儲(chǔ)熱功率
（三）雷諾數(shù)對(duì)堆積床流動(dòng)換熱的影響
基于直徑比為5的堆積床相變儲(chǔ)熱系統(tǒng)，保持相變材料的史蒂芬數(shù)Ste=0.31，當(dāng)換熱流體的進(jìn)口流量分別取220 kg·h?1、260 kg·h?1和300 kg·h?1時(shí)，堆積床的雷諾數(shù)Rep分別為4229.69、4998.72和5767.75。作出不同雷諾數(shù)下堆積床內(nèi)相對(duì)流速的分布，發(fā)現(xiàn)三者的分布幾乎重合，如圖14（a）所示。經(jīng)過(guò)計(jì)算得知，三種雷諾數(shù)下堆積床內(nèi)近壁面區(qū)域換熱流體的最大相對(duì)速度分別是3.40、3.49和3.55，而中心區(qū)域相對(duì)速度的標(biāo)準(zhǔn)差分別是0.40、0.39和0.39，這意味著雷諾數(shù)的增大對(duì)相對(duì)流速的分布幾乎沒(méi)有影響。但由于進(jìn)口速度增大，堆積床內(nèi)的徑向流速發(fā)生成比例的提升，這將增強(qiáng)堆積床內(nèi)的對(duì)流換熱作用，從而加快相變膠囊的儲(chǔ)熱進(jìn)程[圖14 （b）]。
圖14（c）計(jì)算了不同雷諾數(shù)下堆積床相變儲(chǔ)熱系統(tǒng)的儲(chǔ)熱量與平均儲(chǔ)熱功率隨時(shí)間的變化，隨著雷諾數(shù)的增大，系統(tǒng)的儲(chǔ)熱量并不會(huì)改變，但是儲(chǔ)熱時(shí)間縮短了，儲(chǔ)熱功率相應(yīng)提升。雷諾數(shù)分別為 4229.69、4998.72和5767.75時(shí)堆積床的儲(chǔ)熱時(shí)間分別為 97.53 min、89.73 min和83.63 min。但是，雷諾數(shù)的增大也會(huì)引起換熱流體流動(dòng)耗散的增加，從而增大進(jìn)出口的壓降，三個(gè)堆積床的壓降分別為162.84 Pa、213.34 Pa 和270.06 Pa。
圖14 不同雷諾數(shù)堆積床LHS系統(tǒng)的流動(dòng)與熱性能。（a）徑向相對(duì)流速分布；（b）儲(chǔ)熱量與平均儲(chǔ)熱功率；（c）1200 s時(shí)徑向液相分?jǐn)?shù)分布
（四）史蒂芬數(shù)對(duì)堆積床流動(dòng)換熱的影響
與上一部分類(lèi)似，基于直徑比為5的堆積床相變儲(chǔ)熱系統(tǒng)，保持雷諾數(shù)Re=4998.72，當(dāng)換熱流體的進(jìn)口溫度分別取445 ℃、465 ℃和485 ℃時(shí)，堆積床內(nèi)史蒂芬數(shù)Ste分別為0.19、0.31和0.43。圖15（a）為不同史蒂芬數(shù)下堆積床內(nèi)相對(duì)流速的分布，發(fā)現(xiàn)三者的分布重合，這意味著史蒂芬數(shù)的改變對(duì)堆積床內(nèi)流速的大小及分布均無(wú)影響，因此堆積床進(jìn)出口壓降的變化微弱。換熱流體與相變膠囊之間的對(duì)流換熱作用基本不變，但由于史蒂芬的增大代表著進(jìn)口溫度的增大，這將增大換熱流體對(duì)相變膠囊的換熱溫差，從而增強(qiáng)相變膠囊的對(duì)流換熱。如圖15（b）所示，相變材料的儲(chǔ)熱過(guò)程均隨著史蒂芬數(shù)的增大而增大。
史蒂芬數(shù)的增大意味著相變材料存儲(chǔ)的顯熱增加，所以系統(tǒng)的儲(chǔ)熱量隨之增加[圖15（c）]。當(dāng)史蒂芬數(shù)為 0.19時(shí)系統(tǒng)的儲(chǔ)熱時(shí)間為99.20 min，當(dāng)史蒂芬數(shù)增加至 0.31和0.43時(shí)，儲(chǔ)熱時(shí)間變?yōu)?9.73 min和84.7 min。所以，史蒂芬數(shù)的增大不僅能夠提升儲(chǔ)熱系統(tǒng)的儲(chǔ)熱量，還能縮短系統(tǒng)的儲(chǔ)熱時(shí)間，因此能夠很好地提升堆積床相變儲(chǔ)熱系統(tǒng)的熱性能。
圖15 不同雷諾數(shù)堆積床LHS系統(tǒng)的流動(dòng)與熱性能。（a）徑向相對(duì)流速分布；（b）1200 s時(shí)徑向液相分?jǐn)?shù)分布；（c）儲(chǔ)熱量與平均儲(chǔ)熱功率
四、結(jié)論
本文基于壁面效應(yīng)引起的徑向孔隙率振蕩分布，建立了三維堆積床相變儲(chǔ)熱系統(tǒng)模型，通過(guò)將堆積床中相變膠囊溫度變化的模擬值與實(shí)驗(yàn)對(duì)比，驗(yàn)證了模型的可靠性。本文研究了徑向孔隙率振蕩分布對(duì)換熱流體的徑向流速、相變材料的徑向溫度以及相變材料徑向液相分?jǐn)?shù)的影響。此外，還研究了不同的無(wú)量綱參數(shù)（如直徑比、雷諾數(shù)和史蒂芬數(shù)）對(duì)換熱流體和相變材料的影響，主要結(jié)論如下。
（1）通過(guò)沿堆積床徑向截取了一系列的圓柱截面，直觀地解釋了徑向孔隙率振蕩分布的原因。當(dāng)徑向圓柱截面處的相變膠囊為點(diǎn)接觸時(shí)，形成的空隙會(huì)造成高孔隙率，而當(dāng)圓柱截面貫穿相變膠囊時(shí)，相變材料占據(jù)圓柱面上大部分區(qū)域，從而導(dǎo)致低孔隙率。
（2）徑向孔隙率的振蕩分布導(dǎo)致了換熱流體流速在徑向上的不均勻性，在近壁面區(qū)域，換熱流體的流速急劇變化；而在中心區(qū)域流速的分布更為緩和，與徑向孔隙率的分布一致。這也使得孔隙率大的地方相變材料熔化更快。
（3）隨著堆積床直徑比的增大，換熱流體在近壁面區(qū)域的最大流速增大，但是由于相變膠囊的堆積變得更隨機(jī)和致密，換熱流體在中心區(qū)域的流速能夠變得更為均勻。另外，直徑比的增大能夠縮短系統(tǒng)的儲(chǔ)熱時(shí)間并增大平均儲(chǔ)熱功率，但進(jìn)出口壓降會(huì)急劇增大。
（4）換熱流體雷諾數(shù)的增大能夠同比例地增大流速，從而加速相變材料的儲(chǔ)熱過(guò)程，并增大系統(tǒng)的平均儲(chǔ)熱功率。但是進(jìn)出口壓降的大幅提升表明，通過(guò)增大雷諾數(shù)來(lái)提升系統(tǒng)熱性能需在泵功允許的范圍內(nèi)。
（5）相變材料史蒂芬數(shù)的增大對(duì)徑向流速的大小和分布基本無(wú)影響，但是換熱流體與相變材料之間溫差的增加強(qiáng)化了傳熱過(guò)程，并能夠提升系統(tǒng)的總儲(chǔ)熱量和平均儲(chǔ)熱功率。
注：本文內(nèi)容呈現(xiàn)略有調(diào)整，若需可查看原文。
改編原文：
H. B. Liu, C. Y. Zhao.Effect of Radial Porosity Oscillation on the Thermal Performance of Packed Bed Latent Heat Storage[J].Engineering,2021,7(4):515-525.
