多孔介質結構對儲層內流動和換熱特性的影響

汪健生，張輝鵬，2，劉雪玲，2，傅煜郭，2，朱劍嘯，2

（1 中低溫熱能高效利用教育部重點實驗室，天津大學，天津 300350；2 天津大學機械工程學院地熱研究培訓中心，天津 300350）

近年來全球能源格局轉型加快，清潔能源利用比例不斷提高，為解決可再生能源系統在時間和空間上的不匹配特性[1]，儲能技術越來越受到重視。地下儲層受大氣環境波動的影響小，溫度恒定，是一種優良的儲能場所[2-3]。由于含水層儲能的存儲容量高且成本低，利用含水層進行季節性的儲能已在國內外受到廣泛關注。對含水層儲能的研究主要集中于宏觀尺度下儲能系統的熱提取性能，研究者利用均質化模型獲得儲層的平均換熱和流動參數。Ganguly等[4]模擬計算了均質含水層內的瞬態溫度分布，分析了不同流量下的熱提取量。Liu等[5]研究了回灌水溫和孔隙度對儲層內換熱特性的影響，發現孔隙率的影響很大。Ghaebi 等[6]針對一含水層住宅供暖項目，研究了儲能系統在不同運行方案下的熱回收系數。Jiang等[7]改進了抽灌井的過濾結構，研究發現半過濾井結構能夠減緩熱突破，提高傳熱系統的性能系數。但由于含水層儲能對地下環境要求高且可能造成地下水污染，其利用和發展受到限制。通過構建人工儲層，可降低對選址的要求，并減小對地下環境的污染。人工儲層內的填充結構是影響系統性能的關鍵因素。

目前對人工構建地下多孔介質儲層的儲能研究和利用尚不多見。但是在化工工業中，研究人員通過優化填料床內的顆粒填充形狀和填充方式提高換熱效率，加快反應。Calis等[8]模擬了不同粒徑催化劑顆粒填充層內空氣的流動特性，發現管徑和顆粒直徑比低的填料床壓降更小。Alkhalaf 等[9]研究發現顆粒間存在間隙或顆粒重疊時，由于顆粒間的接觸方式不同，影響填料床內流體的交叉混合，進而影響傳熱效果。Gunjal 等[10]比較了4 種不同堆積方式的周期性球體陣列中的流動和傳熱特性，發現了顆粒排列形式的影響較大。Qi等[11]針對非等溫流體的模擬研究，發現填料床的孔隙率是影響流體流動和傳熱的重要因素。考慮到隨機填料床內流體的壓降遠高于結構填料床，Yang等[12]對比研究了均勻和非均勻填料床內的流動和換熱特性，發現選擇合適的填料形式和顆粒形狀，可以顯著降低壓降，提高結構填料床的整體傳熱性能。Fu等[13]實驗研究了非均質多孔介質內的流動傳熱特性，結果表明局部孔隙的差異是造成流體滲流和傳熱差異的主要原因。Halkarni 等[14-15]通過實驗得到了不同粒徑床層的體積對流傳熱系數，發現采用粒徑小的顆粒填充，可提高多孔介質內的體積傳熱系數。Aksornkitti 等[16]研究了熱水在二維顆粒填料床入滲過程中的熱濕傳遞特性，發現流量和粒徑對溫度的影響較大。以上的研究主要集中在球形顆粒填充的多孔介質內的流動和傳熱性能的研究，表明填充顆粒的粒徑、顆粒間間隙以及堆積方式，對多孔介質內的流動和傳熱都有較強的影響。對顆粒的幾何形狀發生變化后，填充多孔介質內的流動和傳熱性能的研究較少。通過構建人工儲層進行地下儲能時，采用球形顆粒填充的儲層蓄水能力較小。為了綜合考慮儲層的傳熱特性及蓄水能力，有必要對非結構化的顆粒填充的多孔介質內的流動和傳熱特性進行探究。為此，構建了非均勻顆粒、十二面體梯度開孔和二十面體梯度開孔的3 種填充顆粒結構，采用共軛傳熱模型，對填充不同顆粒結構的多孔介質通道代表性體積單元內的流動和換熱進行了數值模擬，獲得不同填充結構下的對流傳熱系數和滲流性能，并對其流動和傳熱特性進行綜合評價。通過研究可為人工儲層的構建提供適合的填充物結構，為人工儲層內流體和固體骨架間的對流換熱表面傳熱系數的確定和選取提供可參考的數據。

1 數值模擬過程

1.1 物理模型

地下含水層儲能通常采用雙井系統，回灌井與生產井周圍的儲層結構直接影響含水層的儲能效率和熱提取效果。本研究基于界面表面積密度和孔隙率差異分別構建了3 種多孔介質，即非均勻顆粒、十二面體梯度開孔和二十面體梯度開孔（圖1）。其中圖1(a)填充的非均勻顆粒包括球形以及橢球形，顆粒半徑分布在0.5～0.6mm 之間，圖1(b)和(c)中分別填充邊長均為1mm 的十二面體與二十面體開孔結構，沿著滲流方向孔徑遞減。3種多孔介質填充通道的幾何尺寸見表1。多孔介質骨架為巖石，其熱物性參數見表2。

表1 3種多孔介質填充通道的幾何尺寸

表2 多孔介質骨架材料的物性參數

圖1 不同多孔介質儲層的物理模型

1.2 控制方程和計算方法

首先，基于研究做出以下假設：①流體為不可壓縮流體；②儲層內部不發生生物化學反應；③固體的熱導率為各向同性；④在整個傳熱過程中不考慮輻射換熱。

基于上述假設流體在多孔介質儲層內的流動控制方程如式(1)和式(2)所示。

式中，ρf、μ、p分別為流體的密度、動力黏度和壓力；u為流體速度矢量；I為單位張量；F為體積力。

固體域和流體域的傳熱控制方程如式(3)和式(4)所示。

式中，ks、kf分別為固體和流體的熱導率；cp,f為流體定壓比熱容；Q為流體與固體間的耦合換熱源。水作為儲層內的換熱介質，其物性參數可描述為溫度的函數［式(5)～式(8)］[17]。

研究中使用共軛傳熱模型對不同多孔介質儲層內的對流換熱進行直接模擬，非等溫流體流動與固體和流體傳熱相耦合，通過有限元法對控制方程進行離散化處理，并使用向后微分公式（BDF）求解器對時間變量求解，它是一種向后差分的隱式求解算法，根據物理場對計算收斂容差進行優化控制。

1.3 初始條件和邊界條件

設置的邊界條件和初始條件如下。

（1）儲層入口：Tf=Tin，u=uin

鑒于以上的三種原因，我們特制定如下的改進措施，加大畜牧局檔案管理工作的力度，把這項工作力爭做到實處，取得應有的效果。

（2）儲層出口：p=pout

（3）四周壁面：T=Tw，-n·ρfu=0

（4）骨架壁面：-n·ρfu=0

（5）初始狀態：u=u0，p=p0，Ts=Tw，Tf=T0

1.4 網格無關性檢測和模型驗證

采用四面體網格單元劃分計算域，在填充結構附近對網格進行了加密，表3給出了非均勻顆粒結構在不同網格數下的流體出口溫度，使用Richardson外推法[18]進行了網格獨立性驗證，R小于1 時，數值結果被認為是收斂的。采用同樣的方法對其他兩種結構網格進行了檢驗，最終所使用的網格數見表4。

表3 非均勻顆粒多孔介質在不同網格數下的流體出口溫度（Re=1, Tw=290.65K）

表4 不同多孔介質模型的網格數

為了驗證目前計算模型和方法的可靠性，與Romkes 等[19]所研究的類似問題進行了對比。所采用的物理模型如圖2 所示，8 個均勻的球形顆粒有序堆積在矩形通道內，通道高度與顆粒直徑之比為1，通道壁面絕熱，顆粒壁溫度恒定，空氣為冷卻流體，其進口溫度和速度恒定。將模型求解得出的顆粒壁面與流體間的平均努塞爾數與文獻中的數值模擬和實驗結果進行對比。如圖3所示，計算結果與Romkes 等的實驗以及模擬結果都具有相同的變化趨勢，在低Re（Re<25）下計算結果與Romkes等的結果吻合較好，如Re為10時，與Romkes的計算結果偏差為2.6%。隨著Re的增大，慣性力的影響逐漸增強，多孔介質內的流動將由黏性占主導的Darcy流轉變為Forhheimer流，本文采用的Darcy模型計算偏差會增大。但本文研究涉及的Re小于10，采用該計算模型和方法進行模擬是可靠的。

圖2 Romkes等所使用的物理模型

圖3 計算結果與Romkes等實驗以及模擬結果的對比

2 結果和討論

2.1 流動和換熱關系式

通過對填充單元流動和換熱的模擬計算，得到了不同多孔介質結構的單位壓降（?p/?x）、摩擦系數（f）、骨架壁面與流體間的平均對流傳熱系數（hsf）、平均努塞爾數（Nusf）和總體換熱效率（η），其計算如式(9)～式(17)。

2.2 多孔介質儲層內的流動特性

本節主要討論不同的多孔介質骨架結構儲層內的流動特性。非均勻顆粒、十二面體梯度開孔和二十面體梯度開孔多孔介質儲層內的截面速度分布以及空間流線分布分別如圖4(a)～(c)所示。流體的入口溫度為278.15K。從空間流線分布可以看出，在多孔介質通道內，流體在骨架孔隙內的混合和擾動增強，進而影響流體與固體之間的換熱和壓降。在相同Re和Tw下，模擬結果表明非均勻顆粒通道內的體積平均流速最大，十二面體梯度開孔通道內的體積平均流速最小。非均勻顆粒多孔介質通道內x-y截面速度分布和空間流線分布以及骨架孔隙處的局部放大如圖4(a)所示，x-y截面流速波動明顯，較大的流速主要集中在孔隙中心區域，由于壁面黏性效應，骨架壁面附近的流動阻力較大，局部流體速度較小，遠離壁面后黏性作用減小速度遞增；另外，流體在多孔介質通道內繞掠顆粒，掠過前半部時，流動截面縮小，流速增加，壓力遞降，而在后半部流動截面增加，壓力回升，產生了與滲流方向相反的回流；同時從局部孔隙放大圖可以發現流體在孔隙區域的混合和擾動較強進一步增大了流速。在圖4(b)中，十二面體梯度開孔多孔介質通道內較大的流速主要集中于開孔區域附近而中心區域流速較小，這主要與十二面體開孔結構的開孔位置有關，觀察流線分布可以看出流體進入通道后優先選擇最近的開孔區域繞流，由于流動截面的減小在開孔位置處的流速急劇增加，當進入骨架中心區域后流動截面增大，流速下降；另一方面由于開孔數較少，從y-z截面的速度和流線分布可以看出，中心區域的混合和擾動較弱。與圖4(b)相比，圖4(c)中二十面體梯度開孔結構中心區域的速度較大，改變開孔位置角度以及增加開孔數后流體在中心區域的混合和擾動增強，多方向進行繞流，在中心區域產生了回流和漩渦，這有利于強化流體與多孔骨架之間的傳熱。

圖4 不同多孔介質結構下的截面速度分布和空間流線分布（Re=1, Tw=318.15K）

不同多孔介質結構以及壁面溫度Tw下流體單位壓降的變化如圖5所示，Re范圍在0.01～1，Tw為290.65～318.15K。從圖中可以看出，非均勻顆粒多孔介質通道內流體的?p/?x最大且隨著Re的增大增加幅度也最大，而十二面體梯度開孔多孔介質通道內的?p/?x最小且變化幅度最小，這與局部流速變化相吻合。Re相同時，十二面體梯度開孔結構空間流線分布相對均勻，因此壓降較低，非均勻顆粒結構由于外部繞流彎曲度最大，所以流動阻力更大，壓降最高。隨著Tw的升高，不同多孔介質結構的?p/?x都有所下降，由于低Re下黏性力占主導，溫度升高后流體的黏度下降。圖6給出了不同多孔介質結構以及壁面溫度下摩擦系數的變化。隨著Re的增加，f先是急劇減小然后緩慢下降。Tw相同時，非均勻顆粒結構的f最大，十二面體梯度開孔結構的f最小，這與圖5 中所示的壓降變化相符。隨著Tw的增大，f逐漸減小，在模擬的Re范圍內二十面體梯度開孔與非均勻顆粒兩種多孔介質的f較大且相差較小，而十二面體開孔多孔介質的f則明顯較小，這與多孔介質的結構特性有關，流體在非均勻顆粒多孔介質通道內的迂曲度最大而在十二面體開孔結構通道內的迂曲度最小。以上分析能夠得出，多孔介質結構對流動特性影響顯著，在研究的3種多孔介質結構中，十二面體梯度開孔結構的流動阻力最小。實際應用中選擇合適的多孔介質結構可以降低流體的單位壓降，從而減小流動阻力和流動功耗。

圖5 不同結構和壁面溫度下單位壓降隨Re的變化

圖6 不同結構和壁面溫度下摩擦系數隨Re的變化

2.3 多孔介質儲層內的換熱特性

在多孔介質骨架與通道壁面溫度高于流體溫度時，分析了流體流過多孔結構通道的換熱特性。在Re=1、Tw=318.15K 的工況下，當換熱達到穩定后，不同多孔介質結構的x-y中心截面和距離入口x=0.55mm處y-z截面的溫度分布如圖7所示。比較圖7(a)～(c)的x-y截面云圖可以看出沿著滲流方向不同多孔介質結構整體的溫度分布相似，均呈現出壁面附近溫度較高而中心區域溫度較低的分布。這是因為通道壁面溫度恒定，由于多孔介質骨架與通道壁面之間的導熱作用，通道壁面附近固體骨架的溫度也較高，流體與通道壁面和多孔介質骨架之間的換熱較強。而在通道中心區域由于流體與固體骨架之間的換熱作用骨架溫度降低較快，固體骨架與流體的溫度均較低。對比圖7(a)～(c)的y-z截面溫度云圖可以看出，距入口相同距離處二十面體開孔結構和非均勻顆粒結構中心區域的流體溫度較高，而十二面體開孔結構中心區域溫度則較低，說明十二面體開孔結構流體的換熱效果最差。這與之前的流動分析結果相吻合，十二面體梯度開孔結構在中心區域的混合和擾動較弱，流速較小，換熱較弱；而其他兩種多孔結構通道內流體繞掠骨架過程中在中心區域的混合和擾動較強，因此流速較大，換熱更強。

圖7 不同多孔介質結構下x-y中心截面和距入口0.55mm處y-z截面的溫度分布（Re=1，Tw=318.15K）

不同多孔介質內流體與骨架間的平均對流傳熱系數的變化如圖8所示。可以看出非均勻顆粒的hsf高于其他兩種多孔介質，在Re增大時更加明顯，十二面體梯度開孔的hsf最小，換熱效果最差。Re一定時增大Tw，不同多孔介質的hsf均得到提高，Re越大，升高壁面溫度后對流傳熱系數的增大越明顯。這是由于一方面Tw升高增大了流體與壁面間的換熱溫差，另一方面流體物性隨溫度發生變化，溫度升高后流體的熱導率增大，黏度和比熱容減小，影響流體與壁面之間的換熱。圖9給出了不同多孔介質內Nusf的變化，其中非均勻顆粒結構的Nusf值最大，對流傳熱性能優于另外兩種多孔介質。在Re=0.6～0.8之間，十二面體梯度開孔結構和二十面體梯度開孔結構的Nu值存在交叉，在交叉點之前即Re較小時，二十面體開孔結構的換熱能力占優，隨著Re的增大，二十面體梯度開孔多孔介質的Nusf增加幅度逐漸減緩，而十二面體梯度開孔多孔介質的Nusf則一直保持較高的上升幅度，因此換熱性能逐漸優于二十面體梯度開孔結構。當壁面溫度變化時，由于溫度的改變對流體物性的影響，交叉點也隨之發生變化。

圖8 不同多孔介質與壁面溫度下hsf隨Re的變化

圖9 不同多孔介質與壁面溫度下Nusf隨Re的變化

2.4 多孔介質儲層的總換熱效率

在多孔介質傳熱中，增大流體與骨架之間的換熱面積（即固體骨架的比表面積），可以提高換熱效果，但是流體在多孔介質內部的流動阻力也會增大，造成系統運行過程中循環泵電耗增大。文獻中大多未考慮多孔介質內的流動阻力，僅對地下儲層的換熱性能進行評價[4-6]，不能對儲能系統進行全面評價。為了綜合評價流體在地下儲層流動過程中的換熱效果和流動阻力，本文使用總換熱效率η，即單位壓降的換熱量來評價儲層內的整體換熱性能。不同填充結構和壁面溫度下總換熱效率η的變化如圖10 所示，可以看出隨著Re的增加不同多孔介質儲層的η均略有降低。相同壁面換熱條件下，十二面體梯度開孔多孔介質的η遠高于其他兩種結構，非均勻顆粒與二十面體梯度開孔的η相對較低且比較接近。這表明非均勻顆粒多孔介質和二十面體梯度開孔多孔介質在強化傳熱的同時流動阻力增加較大。所以在3種多孔介質中，盡管十二面體梯度開孔多孔介質的換熱性能低于其他兩種結構，但由于摩擦系數和單位壓降遠低于其他兩種結構，綜合考慮換熱和流動阻力，其總換熱效率最大，綜合性能最好。當Re一定時，隨著Tw的增大，不同多孔介質的η均增大，其中十二面體梯度開孔的增加幅度最大。增大Tw后，一方面冷流體與高溫骨架之間的換熱溫差增大，換熱能力增強，另一方面由于流體黏度隨溫度升高而減小，流動阻力減小，因此不同多孔介質儲層的總傳熱效率隨著壁面溫度的增大均增大。

圖10 不同多孔介質結構和壁面溫度下η隨Re的變化

3 結論

對非均勻顆粒、十二面體梯度開孔和二十面體梯度開孔3種多孔介質內的流動和換熱特性進行了模擬研究，對比分析了3種結構內流體的單位壓降、摩擦系數、對流傳熱系數、平均努塞爾數及綜合換熱效率，并探討了Re以及壁面溫度Tw的影響。研究得到了3 種結構填充的多孔介質內的流動特性、傳熱特性以及總換熱效率的大小，在人工構建地下含水儲層時，可指導選取合適的填充材料的幾何結構。主要結論如下。

（1）在相同的Re和Tw下，3種多孔介質中十二面體梯度開孔的?p/?x和f最小，非均勻顆粒的?p/?x和f最大。減小Re和增大Tw均可減小單位壓降，從而減小流動阻力及循環泵功耗。

（2）在相同的Re和Tw下，3種多孔介質中非均勻顆粒的hsf與Nusf最高，十二面體梯度開孔的hsf最小。二十面體梯度開孔多孔介質與十二面體梯度開孔多孔介質的Nu隨Re的變化存在交叉，Re較小時二十面體梯度開孔的Nusf較大，Re較大時十二面體梯度開孔的Nusf較大。交叉點出現在Re為0.6～0.8之間，交叉點的Re會受Tw的影響而變化。

（3）綜合考慮壓降和換熱性能的評價指標總換熱效率η，十二面體梯度開孔多孔介質的最高，二十面體梯度開孔與非均勻顆粒的比較接近。隨著Re的增大不同結構的η均有所下降，Tw升高后不同結構的η得到明顯改善。