中深層地埋管換熱器取熱穩(wěn)定性及熱影響半徑

黃 帥，朱 科，董建鍇，李 驥，姜益強，方肇洪

(1.哈爾濱工業(yè)大學 建筑學院，哈爾濱 150006；2.寒地城鄉(xiāng)人居環(huán)境科學與技術(shù)工業(yè)和信息化部重點實驗室(哈爾濱工業(yè)大學)，哈爾濱 150090；3.山東建筑大學 熱能工程學院，濟南 250101；4.中國建筑科學研究院有限公司，北京 100101)

地熱能資源豐富，具有分布廣、清潔環(huán)保、儲量大、穩(wěn)定可靠等特點[1]。在建筑行業(yè)地熱能開發(fā)利用過程中，淺層地源熱泵系統(tǒng)應用較為廣泛[2]，但其對空調(diào)全年冷熱負荷要求較高，長期運行后系統(tǒng)性能會顯著下降[3]，并且淺層地埋管換熱器埋管占地面積較大。近年來，中深層地源熱泵系統(tǒng)作為新興的地熱能利用形式，通常其埋管深度在1 000～3 000 m，具有地埋管換熱器占地面積小、對全年冷熱負荷平衡要求低、系統(tǒng)取熱性能穩(wěn)定等優(yōu)勢[4]，受到了廣泛關(guān)注。

與傳統(tǒng)淺層地埋管換熱器相比，中深層地埋管供熱技術(shù)只在近年來得以迅速發(fā)展[5]，并且其發(fā)展現(xiàn)狀是工程實踐先于理論研究，有很多問題尚待于研究[6]。其中，在中深層地埋管換熱器取熱特性方面，Li等[7]對影響中深層地埋管換熱器取熱特性的因素進行分析，研究結(jié)果表明增大外管管徑和流速有利于提升中深層地埋管換熱器的取熱量；Liu等[8]對中深層地埋管換熱器管內(nèi)循環(huán)工質(zhì)的流速進行了優(yōu)化，并提出最佳流速為0.7 m/s；Fang等[9]發(fā)現(xiàn)增大中深層地埋管換熱器的內(nèi)套管熱阻能夠提升其換熱能力。而對中深層地埋管換熱器周圍的巖土熱恢復特性等問題研究較少，其中，賈林瑞等[10]研究了單取熱工況地下巖土溫度熱響應，發(fā)現(xiàn)1 400～1 800 m處的巖土溫度變化率最快；Cai等[11]對間歇運行模式下中深層地埋管換熱器取熱穩(wěn)定性進行了研究，發(fā)現(xiàn)不同運停比運行時的巖土溫度熱恢復效果較好；Song等[12]研究結(jié)果表明巖土經(jīng)過為期8個月的熱恢復期，其溫度基本能恢復至初始狀態(tài)。但上述研究主要是針對單個供熱周期，其目的是量化設(shè)計參數(shù)對中深層套管式地埋管換熱器取熱性能的影響以及揭示巖土的熱恢復規(guī)律。考慮到中深層地熱地埋管供熱系統(tǒng)需長期運行，因此，為明確利用該系統(tǒng)供熱的可持續(xù)性，需要對中深層套管式地埋管換熱器長期取熱的穩(wěn)定性進行相關(guān)研究。

此外，中深層地埋管換熱器取熱對周圍巖土的熱擾與諸多參數(shù)相關(guān)，如運行時間、管徑、循環(huán)流量、回填材料導熱系數(shù)及巖土熱物性等。在工程實踐過程中，為解決中深層地埋管換熱器鉆孔間距不當造成的熱短路以及土地資源浪費等問題，首先要確定中深層地埋管換熱器取熱對周圍巖土的熱擾與上述哪些參數(shù)相關(guān)性較大，因此，需要對上述參數(shù)做敏感性分析。

基于此，本文建立了中深層地埋管換熱器數(shù)值傳熱模型，基于有限差分法將控制方程離散求解，并利用項目實測數(shù)據(jù)進行驗證，在此基礎(chǔ)上研究長期運行下中深層套管式地埋管換熱器的出水溫度以及巖土的溫度變化。此外，對影響中深層套管式地埋管換熱器熱擾的運行參數(shù)進行了敏感性分析，探尋中深層地埋管換熱器的熱影響半徑與相關(guān)參數(shù)的表達式，為日后的工程實踐提供理論指導。

1 數(shù)值模型建立

1.1 物理模型

中深層同軸套管式地埋管換熱器的取熱原理如圖1所示，在循環(huán)水泵的驅(qū)動下循環(huán)工質(zhì)(水)由套管外部環(huán)腔注入，經(jīng)地埋管換熱器周圍高溫巖土層加熱后通過內(nèi)套管反向流動，流出地面后進入熱泵機組。

圖1 中深層地埋管換熱器的換熱原理Fig.1 Heat exchange principle of MDBHE

中深層地埋管換熱器的取熱過程主要包括內(nèi)套管流體與外部環(huán)腔流體的對流換熱、環(huán)腔流體與外管壁的對流換熱以及外套管與周圍巖土的導熱，其換熱過程較為復雜，并且地下巖土環(huán)境具有不確定性，因此做出以下假設(shè)：

1)忽略地下滲流作用的影響，將巖土中傳熱視為單純的導熱問題[13]；

2)忽略地表溫度波動的影響[14]；

3)假定數(shù)值模擬區(qū)域的徑向邊界處的溫度分布不受地埋管換熱器的影響；

4)假定同軸套管式換熱器內(nèi)循環(huán)體主要以對流傳熱，忽略其軸向的熱傳導。

1.2 控制方程

基于以上假設(shè)，每層巖土的導熱方程可寫為[15]

(1)

式中：ak為不同巖土層的熱擴散率，m2/s；Ts為巖土層溫度，℃；Δτ為時間步長，s；r、z分別為徑向長度與豎直方向長度，m。

當管內(nèi)流體工質(zhì)以圖1所示的方式流動時，外部環(huán)腔內(nèi)流體的能量方程為

(2)

內(nèi)管流體的能量方程為

(3)

式(2)、(3)描述的是在中深層地埋管換熱器正常取熱階段下的能量方程，當系統(tǒng)停運時，外部環(huán)腔內(nèi)流體及內(nèi)管流體的流速為零，此時地埋管周圍巖土處于熱恢復狀態(tài)。與正常取熱階段下的能量方程相比，熱恢復階段能量方程中的系數(shù)C為零。

其中：

式中：C=MCf為循環(huán)液的熱容流量,kJ/(s·K)；db為鉆孔半徑,m；di、do分別為套管的內(nèi)徑和外徑，m；ρc、ρ1c1、ρ2c2、ρgcg分別為循環(huán)水比定容熱容，外管比定容熱容、內(nèi)管比定容熱容、回填材料比定容熱容，J/(m3·K)。此外，外管內(nèi)循環(huán)水和鉆孔壁之間的熱阻R1以及內(nèi)外管循環(huán)水之間的熱阻R2由下式求得：

(4)

(5)

式中：kg、kp1、kp2分別為回填材料、外管和內(nèi)管的導熱系數(shù)，W/(m·K)；h1、h2分別為套管外管和內(nèi)管的對流換熱系數(shù)，W/(m2·K)。

1.3 初始及邊界條件

在任一深度的地層的初始地溫可以表達為下式[9]：

(6)

式中：Hj為第j層地層底部的坐標；qg為大地熱流，W/m2；ta為地表溫度，℃；ha為空氣與地表的對流換熱系數(shù)，W/(m2·K)；km為巖土的導熱系數(shù)，W/(m·K)。

外部環(huán)腔流體與內(nèi)管流體能量方程的邊界體條件為：

(7)

Tf1=Tf2,Z=H

(8)

式中：Q為取熱功率，kW；C=MCf為循環(huán)液的熱容流量,kJ/(s·K)；Z為埋管深度，m。

2 數(shù)值模型求解及驗證

2.1 模型求解

在非穩(wěn)態(tài)數(shù)值傳熱模擬中，時間步長Δτ、縱向步長Δz及徑向步長對模擬計算速度及結(jié)果的精確度均有一定影響。為了保證結(jié)果的準確性需要對其做無關(guān)性驗證，由圖2可知Δτ可取900 s，此外，考慮徑向較遠處熱流趨近于零，鉆孔附近的熱流及溫度梯度較大，故以鉆孔壁為界采用變步長，即徑向步長按照等比級數(shù)變化(等比系數(shù)為1.2)，徑向離鉆孔第1個節(jié)點距離為0.168 m,徑向共取40個節(jié)點，縱向步長Δz取10 m。由圖3可知節(jié)點數(shù)量可取7 700個，滿足計算精度要求。其余模擬所需的基礎(chǔ)物理參數(shù)見表1，巖土熱物性參數(shù)除特別說明外取下文中的案例1。

表1 模擬所需的物理參數(shù)Tab.1 Physical parameters required for simulation

圖2 時間步長無關(guān)性驗證Fig.2 Verification of time step independence

圖3 節(jié)點數(shù)量無關(guān)性驗證Fig.3 Verification of node number independence

在中深層地埋管換熱器與周圍巖土層的耦合傳熱問題中，可將其看作圓柱體的二維瞬態(tài)導熱。對于二維瞬態(tài)耦合傳熱問題，若采用隱式格式差分方法得到非邊界處節(jié)點方程中的未知數(shù)通常有5個，采用迭代法求解耗時較長[16]。因此，本文基于微元體熱平衡和交替時間步長法[17]建立節(jié)點方程，離散網(wǎng)格如圖4所示。

圖4 模型離散網(wǎng)格示意Fig.4 Schematic of discretization mesh for numerical model

2.2 模型驗證

為了驗證建立的數(shù)學模型及所編程序的準確性，作者所在團隊于2019年11月至2020年1月期間對西安市某住宅工程進行了現(xiàn)場實測，測試參數(shù)主要包括中深層地埋管換熱器進出水溫度、機組及循環(huán)水泵耗電量、用戶側(cè)進出水溫度及循環(huán)水泵耗電量，測試期間機組的平均COP達4.7，限于篇幅，其余內(nèi)容詳見文獻[7]。其中，測試期間單孔地熱井平均換熱量為270 kW，循環(huán)水平均流量為25.5 m3/h，中深層地埋管換熱器進出水溫差約為9.2 ℃，將該地熱井運行參數(shù)導入所編制的程序，進而得出模擬地源側(cè)進出口溫度，選取650 h后的測試數(shù)據(jù)與模擬值進行對比，地源側(cè)進出口溫度的模擬值與實測值對比如圖5所示。從圖5中出水溫度可知，除個別波動點外，最大誤差為6.01%，進而驗證了數(shù)值傳熱模型及所編程序的正確性，為下文相關(guān)分析奠定了基礎(chǔ)。

圖5 地源側(cè)進、出口溫度實際值與模擬值比對Fig.5 Comparison of actual and simulated temperature at inlet and outlet of ground source side

3 數(shù)值模擬結(jié)果與分析

考慮到中深層地熱地埋管供熱系統(tǒng)需長期運行，因此，為明確利用該系統(tǒng)供熱的可持續(xù)性，需要對中深層套管式地埋管換熱器長期取熱的穩(wěn)定性以及對其周圍巖土的熱擾進行相關(guān)研究。

3.1 中深層套管式地埋管換熱器的穩(wěn)定性

3.1.1 長期出水溫度

對中深層套管式地埋管換熱器的長期出水溫度進行為期15 a的模擬計算，此時的延米換熱量為100 W/m，流量為28 m3/h，其余參數(shù)見表1。

圖6所示為中深層套管式地埋管換熱器在15 a運行期間的出水溫度。從圖6中可以看出，在每個供熱周期內(nèi)的中深層套管式地埋管換熱器的出水溫度變化較大，以第1個供暖期為例，初始階段的出水溫度為35.6 ℃，采暖季結(jié)束后的出水溫度為25.9 ℃。從整體看，中深層套管式地埋管換熱器的出水溫度隨著運行年限的增加呈遞減趨勢，并且遞減程度逐漸減小，第1 a、第5 a、第10 a和第15 a采暖季末期的中深層套管式地埋管換熱器的出水溫度分別為25.9、24.4、23.9、23.5 ℃。分別以第1 a、第5 a和第10 a的出水溫度為基準，第5 a、第10 a和第15 a的出水溫度在此基礎(chǔ)上分別下降5.8%、2.0%和1.8%。此外，從圖6中也可以看出，從第11 a開始出水溫度基本處于準穩(wěn)態(tài)階段，第11 a、第12 a、第13 a、第14 a和第15 a采暖季末期的中深層套管式地埋管換熱器的出水溫度分別為23.79、23.72、23.66、23.60、23.54 ℃，表明中深層套管式地埋管換熱器長期運行下的出水溫度較為穩(wěn)定。

圖6 15 a內(nèi)的中深層套管式地埋管換熱器出水溫度Fig.6 The oulet temperature of MDBHE within 15 a

3.1.2 長期巖土溫度變化

對距地埋管換熱器外壁1 m，不同埋管深度的巖土溫度進行為期15 a的模擬計算，此時的延米換熱量為100 W/m，流量為28 m3/h，其余參數(shù)見表1。

從圖7可以看出，在15 a運行期間，距鉆孔壁1 m處不同深度巖土層的溫度呈周期性變化，并且?guī)r土經(jīng)過熱恢復期不能完全恢復到初始狀態(tài)，與初始溫度相比，巖土溫度隨運行年限的增加逐年下降。以z=1 800 m為例進行分析，與巖土初始溫度相比，第1 a溫降比例為2.03%，第5 a溫降比例為4.30%，第10 a溫降比例為5.37%，第13～15 a溫降比例分別為5.78%，5.90%和6.01%，溫降增幅逐年減小。分析原因：在運行參數(shù)不變的條件下，被取熱后的巖土不能完全恢復到初始溫度，會造成下一年地埋管換熱器取熱能力略有下降，進而對其周圍巖土的熱干擾作用減弱，經(jīng)過長時間的周期循環(huán)這種現(xiàn)象會逐漸減弱，最后會趨于相對平衡，此外，經(jīng)過為期15 a運行，最大溫降比例僅為6.50%，說明中深層地埋管換熱器取熱穩(wěn)定性較好，巖土具有較強的熱恢復能力。

圖7 不同埋管深度下徑向1 m處的巖土溫度變化Fig.7 Temperature variation of rock and soil at 1 m in radial direction under different buried depths

3.2 中深層套管式地埋管換熱器的熱擾

Δtr,z,τ=tbnd,z,0-tr*,z,τ

(9)

式中：tbnd,z,0為τ=0時刻徑向邊界溫度，℃；tr,z,τ為τ時刻半徑r處的溫度,℃。

3.2.1 不同運行時間下熱影響半徑

圖8所示為中深層地埋管換熱器運行5、10、15 a后的巖土溫度場以及其不同埋管深度處徑向巖土溫度的變化曲線。以圖8(f)為例進行分析，當運行15 a后，與徑向邊界171.47 m處的巖土溫度相比，當中深層地埋管換熱器的埋管深度分別為1 200、1 400、1 600、1 800、2 000 m時，距離其0.17 m處的巖土溫降分別為20.79、23.07、22.91、23.81、24.51 ℃，相比而言z=1 200 m處的溫降最小，表明中深層地埋管換熱器的取熱區(qū)域主要集中在1 400～2 000 m處的巖土。此外，隨著運行時間的增加中深層地埋管換熱器的熱影響半徑逐漸增大，第5、10、15 a的影響半徑分別為39.87、57.42、82.69 m，以第10 a對應的熱影響半徑為基準，第15 a的熱影響半徑在其基礎(chǔ)上增加44%。分析原因：由前文分析可知經(jīng)取熱后的巖土不能完全恢復到初始狀態(tài)，并且，隨著運行年限的增加，巖土溫度逐漸降低，由于熱傳導的作用，周圍未受影響的巖土會向地埋管中心進行熱量傳遞，故熱影響半徑會隨運行時間的增加而逐漸增大。

圖8 不同運行年限下溫度場分布和不同埋管深度徑向溫度的變化Fig.8 Distribution of temperature field and variation of radial temperature at different buried depths under different operation periods

3.2.2 主要設(shè)計參數(shù)與熱影響半徑的關(guān)系

采用控制變量法分析主要設(shè)計參數(shù)和中深層地埋管換熱器熱影響半徑之間的關(guān)系，其中，主要設(shè)計參數(shù)包括管徑、工質(zhì)循環(huán)流量、回填材料導熱系數(shù)以及延米換熱量。在不同設(shè)計參數(shù)下，系統(tǒng)運行10 a后z=2 000 m處徑向巖土溫度的分布如圖9所示。從圖9可以看出，管徑、工質(zhì)循環(huán)流量以及回填材料導熱系數(shù)的變化對中深層地埋管換熱器周圍巖土層溫度分布影響較小，相對于其他3種設(shè)計參數(shù)，延米換熱量對中深層地埋管換熱器周圍巖土層溫度分布影響較大，隨著延米換熱換熱量的增大周圍巖土層的溫度越低。因此，在進行工程設(shè)計時應采用合理的延米換熱量，文獻[18]對延米換熱量做了詳細研究，并指出延米換熱量不宜大于150 W/m。此外，從圖中可以看出在一定范圍內(nèi)，上述4種設(shè)計參數(shù)雖對中深層地埋管換熱器周圍巖土具有一定的影響，但對其熱半徑的熱影響范圍有限，即中深層地埋管換熱器的熱影響半徑趨于一致。

圖9 不同設(shè)計參數(shù)下運行10 a后徑向巖土溫度的分布Fig.9 Radial temperature distribution of rock and soil under different design parameters after 10 a operation

3.2.3 巖土熱物性參數(shù)和熱影響半徑的關(guān)系

本文主要分析巖土熱物性參數(shù)和熱影響半徑之間的關(guān)系。巖土不同熱物性參數(shù)見表2，此時的延米換熱量為100 W/m，流量為28 m3/h，其余參數(shù)見表1。

在表2所示的巖土熱物性參數(shù)下，圖10所示為z=500 m處巖土的徑向溫度分布。從圖10中可以看出，與徑向邊界171.47 m處的巖土溫度相比，當巖土熱物性參數(shù)分別取案例1～5時，距離中深層地埋管換熱器0.17 m處的巖土溫降分別為10.55、11.61、16.13、14.94、5.58 ℃，由此可知，中深層地埋管換熱器取熱對其周圍巖土溫度場的熱擾主要與巖土本身熱物性參數(shù)相關(guān)。

表2 巖土物理參數(shù)Tab.2 Physical parameters of rock and soil

圖10 不同巖土熱物性參數(shù)下徑向巖土溫度分布Fig.10 Radial temperature distribution along depth of DBHE under different geotechnical thermal physical parameters

案例1～5的巖土單位體積平均熱容分別為2.410×106、1.390×106、2.040×106、2.210×106、2.680×106J/(m3·K)，其中案例2的單位體積平均熱容最小，案例5的單位體積平均熱容最大。當z=2 000 m時，案例2的熱影響半徑最大，而案例5的熱影響半徑最小，其中，案例1～5的巖土平均熱擴散率分別為1.43×10-6、2.01×10-6、8.60×10-7、6.19×10-7、4.13×10-7m2/s，所對應的中深層地埋管換熱器熱影響半徑分別為82.69、99.23、68.91、57.43、27.69 m。因此，中深層地埋管換熱器熱影響半徑與巖土單位體積熱容及熱擴散率等物性參數(shù)相關(guān)性較大，巖土單位體積熱容越小、熱擴散率越大，熱影響半徑也越大。但隨著熱擴散率的增大，中深層地埋管換熱器的半徑增加的速率逐漸減小，當熱擴散率由4.13×10-7m2/s增加到6.19×10-7m2/s時，熱影響半徑由27.69 m增加到57.43 m；而當熱擴散率由1.43×10-6m2/s增加到2.01×10-6m2/s時，熱影響半徑由82.69 m增加到99.23 m，增加率為20%。

為進一步探究中深層地埋管換熱器熱影響半徑與巖土熱擴散率的關(guān)系，圖11所示為兩者擬合關(guān)系式。由圖11(a)可知，兩者并非簡單的線性關(guān)系，具體表現(xiàn)為隨巖土熱擴散率的增加，中深層地埋管換熱器的熱影響半徑在前期增加速率較快，而后期較慢，滿足指數(shù)關(guān)系式的變化規(guī)律。從圖11(b)中可以看出，當熱影響半徑r和熱擴散率a的關(guān)系式滿足r=-151.99×exp[-a/(5.14×10-7)]+98.14時，此時的R2為0.97，能夠較好的表述兩者之間的關(guān)系。

圖11 巖土熱擴散率與中深層地埋管熱影響半徑的關(guān)系Fig.11 Relationship between thermal diffusivity of rock and soil and heat-affected radius of MDBHE

但應當指出：利用中深層地埋管換熱器熱影響半徑確定鉆孔間距并不是最好的方法，熱影響半徑和間距問題并不完全等價，最好的方法是直接求解多根中深層地埋管換熱器之間耦合，而且要求解反問題，即算出取熱量隨間距的變化，當然難度和計算工作量較大。現(xiàn)有文獻的報道多關(guān)于單根中深層地埋管換熱器，但在實際工程中，延米換熱量取值一般要小于150 W/m，單根中深層地埋管換熱器不能滿足現(xiàn)有公共建筑及住宅建筑的供暖需求，往往需要更多根中深層地埋管換熱器來滿足建筑總供熱負荷。因此，為使研究成果更好的服務(wù)實際工程，應研究中深層地埋管換熱器管群，分析管群的布管形式以及地埋管換熱器之間的相互影響，但在現(xiàn)有研究背景下，本研究成果仍有一定的參考意義。

4 結(jié) 論

1)中深層套管式地埋管換熱器長期運行下的出水溫度較為穩(wěn)定。從第11 a開始出水溫度基本處于準穩(wěn)態(tài)階段，第11 a、第12 a、第13 a、第14 a和第15 a采暖季末期的中深層套管式地埋管換熱器的出水溫度分別為23.79、23.72、23.66、23.60、23.54 ℃。

2)中深層套管式地埋管換熱器周圍巖土熱恢復能力較強。對距地埋管換熱器外壁1 m位置處不同埋管深度的巖土溫度進行為期15 a計算，結(jié)果顯示經(jīng)過15 a的運行，與初始巖土溫度相比最大溫降比例僅為6.50%。

3)中深層套管式地埋管換熱器的熱影響半徑隨運行時間的增加而逐漸增大。運行5、10、15 a的影響半徑分別為39.87、57.42、82.69 m，以第10 a對應的熱影響半徑為基準，5 a后的熱影響半徑在其基礎(chǔ)上增加44%。

4)相比于中深層套管式地埋管換熱器外管管徑、回填材料導熱系數(shù)及工質(zhì)循環(huán)流量，延米換熱量對中深層換熱器周圍巖土的溫度分布影響較大。但在一定范圍內(nèi)，上述設(shè)計參數(shù)和中深層套管式地埋管換熱器熱影響半徑之間無明顯關(guān)系。

5)隨巖土熱擴散率的增加，中深層套管式地埋管換熱器的熱影響半徑在前期增加速率較快，而后期較慢，當熱擴散率由4.13×10-7m2/s增加到6.19×10-7m2/s時，熱影響半徑由27.69 m增加到57.43 m；而當熱擴散率由1.43×10-6m2/s增加到2.01×10-6m2/s時，熱影響半徑由82.69 m增加到99.23 m，增長率為20%，滿足指數(shù)關(guān)系式的變化規(guī)律。并且，熱影響半徑r和熱擴散率a的關(guān)系式滿足r=-151.99×exp[-a/(5.14×10-7)]+98.14，其中，R2=0.97。