熱儲分層影響EGS采熱的數(shù)值模擬研究*

黃文博，陳繼良，蔣方明

（中國科學(xué)院廣州能源研究所，中國科學(xué)院可再生能源重點實驗室，先進能源系統(tǒng)實驗室，廣州 510640）

熱儲分層影響EGS采熱的數(shù)值模擬研究*

黃文博，陳繼良，蔣方明?

（中國科學(xué)院廣州能源研究所，中國科學(xué)院可再生能源重點實驗室，先進能源系統(tǒng)實驗室，廣州 510640）

由于巖石構(gòu)造不同、天然裂隙的差異以及壓裂過程的隨機性等因素，增強型地?zé)嵯到y(tǒng)（EGS）人工熱儲通常具有較強的非均質(zhì)性。探究熱儲的非均質(zhì)性對EGS熱開采過程的影響，對EGS性能預(yù)測與分析評價有重要意義。論文考慮到熱儲沿深度方向的非均質(zhì)性，基于等效分層多孔介質(zhì)物理模型，并使用自主開發(fā)的EGS數(shù)值模型，模擬了多個具有分層熱儲EGS的長期運行過程，發(fā)現(xiàn)熱儲深度方向上非均質(zhì)性對熱能的開采影響顯著，而流量分布的不均勻性是導(dǎo)致系統(tǒng)采熱性能下降的主要原因。為了方便分析和評價，我們建立一種新的定量化描述熱儲非均質(zhì)性的方法，然后基于更多的非均質(zhì)熱儲EGS算例結(jié)果，擬合得到EGS采熱性能與熱儲非均質(zhì)性的定量關(guān)系式。

增強型地?zé)嵯到y(tǒng)（EGS）；數(shù)值模擬；等效分層多孔介質(zhì)；非均質(zhì)性

0 前 言

在地下3～10 km的低滲透性巖石中儲存著大量的熱能。為了有效地開發(fā)利用這些熱能，美國拉斯阿莫斯國家實驗室（Los Alamos National Laboratory）在 1970年提出了增強型地?zé)嵯到y(tǒng)（Enhanced Geothermal Systems, EGS）概念——通過水力壓裂、化學(xué)腐蝕等方法增加地下深層巖石的滲透性形成人工熱儲，然后建設(shè)流體循環(huán)系統(tǒng)，經(jīng)由注入井注入冷流體工質(zhì)，其在人工熱儲被加熱后由采出井輸送至地面電廠，發(fā)電后的流體經(jīng)進一步的梯級利用后再灌注到地下熱儲，循環(huán)利用，從而實現(xiàn)深層地?zé)崮艿某掷m(xù)開采[1]。EGS不僅具有常規(guī)水熱型地?zé)峁┠芊€(wěn)定持續(xù)且污染物排放量小的優(yōu)點，同時還具有儲量巨大、分布廣泛的特點，有望成為未來能源中的重要組成部分[1-5]。

人工熱儲的建立是 EGS建設(shè)過程中的關(guān)鍵環(huán)節(jié)，水力壓裂是通常采用的方法。該方法通過往深層巖石中持續(xù)注入高壓液體，壓裂巖石或使巖石中的天然裂縫發(fā)生擴張和延伸，形成結(jié)構(gòu)復(fù)雜的裂隙網(wǎng)絡(luò)，增加巖石的滲透性[1]。由于巖石構(gòu)造不同、天然裂縫的差異以及壓裂過程的各種不確定性因素，熱儲區(qū)域裂隙網(wǎng)絡(luò)往往具有很強的非均質(zhì)性[6-10]。熱儲的裂隙結(jié)構(gòu)直接影響流體在熱儲中滲流換熱過程，是決定熱儲可開采熱能和EGS壽命的關(guān)鍵因素。評估熱儲裂隙網(wǎng)絡(luò)的非均質(zhì)性對EGS采熱性能的影響，可為人工熱儲的壓裂激發(fā)過程提供指導(dǎo)，對提高EGS經(jīng)濟性具有重要意義。

數(shù)值模擬作為一種成本低、效率高并且功能強大的研究方法，在當前的EGS研究工作中得到了廣泛的應(yīng)用。EGS數(shù)值模擬中，根據(jù)對裂隙網(wǎng)絡(luò)處理方式的不同，熱儲模型大致可分為兩類[11,12]：離散裂隙模型（Discrete Fracture Network Model, DFN）和連續(xù)介質(zhì)模型（Equivalent Continuum Model）。離散裂隙模型嚴格地區(qū)分巖石骨架和裂隙，可以很好地體現(xiàn)裂隙開度、裂隙空間分布等因素對流體流動和換熱的影響。但是EGS熱儲尺寸通常為裂隙尺寸的105～106倍，對于實際的EGS熱儲，DFN會因計算量過大難以進行。連續(xù)介質(zhì)模型基于等效體積（Representative Equivalent Volume, REV）平均的假設(shè)，即使用體積平均性質(zhì)來近似描述熱儲局部的各種性質(zhì)，可以在未探明熱儲內(nèi)裂隙網(wǎng)絡(luò)的具體幾何形態(tài)和空間分布的情況下對熱儲進行較為準確的宏觀等效描述，計算量相對較小，是目前EGS數(shù)值模型中主要使用的方法。

在連續(xù)介質(zhì)模型中，根據(jù)對熱儲參數(shù)的設(shè)定方法的不同，等效多孔介質(zhì)又可分為單孔隙度模型（Single-porosity Model）[13-16]、雙孔隙度模型（Double-porosity Model）[17-20]和多孔隙度模型（Multi-porosity Model）[21-26]三種。

由于熱儲位置通常在地下數(shù)千米，獲得精確描述不均勻熱儲的結(jié)構(gòu)信息在目前來看依然十分困難，很多文獻中采用了單孔隙度模型。這種模型將熱儲視為均勻熱儲，認為熱儲中的地?zé)岬刭|(zhì)參數(shù)大致為該區(qū)域的平均值。單孔隙度模型通常用于檢驗EGS數(shù)值模型的適用性[13]，也可以在一定程度上對EGS宏觀運行規(guī)律進行研究[14-16]。由于這種模型忽略了熱儲的非均質(zhì)性，難于真實反映熱儲中的滲流過程，對實際EGS系統(tǒng)運行情況的預(yù)測會有偏差。

為了較好地體現(xiàn)實際熱儲中的結(jié)構(gòu)，一些研究者采用了雙孔隙度模型。雙孔隙度模型仍視熱儲為多孔介質(zhì)，但它將每個熱儲單元劃分為兩類子域：低孔隙率的巖石骨架和高孔隙率的裂隙。在雙孔隙度模型基礎(chǔ)上，多重連續(xù)介質(zhì)（Multiple Interacting Continua, MINC）模型對雙孔隙度模型中區(qū)別明顯的兩個子區(qū)域進一步劃分，避免在局部地區(qū)產(chǎn)生過大的物理量梯度，使得計算更趨穩(wěn)定。以MINC模型為基礎(chǔ)的 TOUGH軟件目前得到了廣泛的應(yīng)用[19,20]，但雙孔隙率模型面臨巖石骨架和裂隙區(qū)域如何劃分、如何根據(jù)野外地質(zhì)勘測和實驗數(shù)據(jù)區(qū)別設(shè)置兩區(qū)域的參數(shù)等困難。

為了更好地把實地勘測的數(shù)據(jù)反映到熱儲模型中，一些研究者采用了多孔隙度模型。Kalinina等[21]、 Shaik等[22]將裂隙方向、間距、開度和長度等重要裂隙參數(shù)都考慮進模型中，進而構(gòu)建出孔隙度、滲透率的空間分布；Vogt等[23,24]嘗試通過 EnKF（Ensemble Kalman Filter）和蒙特卡洛方法（Monte Carlo approach），依據(jù)示蹤劑循環(huán)實驗等實測結(jié)果重構(gòu)出滲透率在熱儲中的分布情況。但是，由于難以獲得足夠多的精確地質(zhì)數(shù)據(jù)，目前多孔隙度模型仍不足以還原真實EGS熱儲復(fù)雜的裂隙結(jié)構(gòu)。

Fourar[25]提出了用等效 分 層 多孔介質(zhì)（Equivalent Stratified Porous Medium）來描述不均勻熱儲。這種模型將熱儲假設(shè)為分層結(jié)構(gòu)，各個流動層具有不同的滲透性，并且認為熱儲中的流動平行于流動層。Radilla等[26]嘗試將該模型用于解釋法國蘇爾士（Soultz）地區(qū)EGS熱儲結(jié)構(gòu)，其計算結(jié)果與示蹤粒子實驗結(jié)果吻合很好。Bl?cher等[27]將分層多孔巖石儲層模型應(yīng)用于模擬深層地?zé)醿Τ貎?nèi)水熱過程。最近，Luo等[28,29]還將分層熱儲模型用于CO2-EGS采熱過程[28]和CO2地質(zhì)存儲[29]的模擬。

本文將EGS熱儲考慮成多個單孔隙度流動層，通過設(shè)定各流動層的孔隙度、滲透率等宏觀參數(shù)來模擬構(gòu)建沿深度方向的不均勻熱儲，基于對EGS的長期采熱過程的數(shù)值模擬結(jié)果，研究熱儲深度方向的非均質(zhì)性對EGS采熱性能的影響。本文將在詳細分析熱儲滲流分布對EGS采熱性能影響的基礎(chǔ)上，定義合理的描述熱儲非均質(zhì)性參數(shù)；同時，將提出相對采熱率的概念并采用其評價不均勻熱儲EGS的采熱性能；另外，還將采用隨機方法構(gòu)建30個分層不均勻熱儲并模擬EGS長期運行，繼而通過分析數(shù)值結(jié)果確定分層熱儲結(jié)構(gòu)下EGS采熱性能與熱儲非均質(zhì)性的定量關(guān)系。

1 數(shù)學(xué)-物理模型

1.1 物理模型

考慮典型EGS雙井布局，EGS的地下部分由井、人工熱儲和熱儲周圍的干熱巖三個部分組成。在數(shù)值計算中，忽略熱儲周圍巖石的滲透性，注入井和生產(chǎn)井均視為開放式流體通道。如圖1所示，Region 1代表分層熱儲，Region 2為熱儲周圍不可滲透的巖石，Region 3為注入井和生產(chǎn)井。

該模型對熱儲進行了以下假設(shè)：（1）熱儲視為分層結(jié)構(gòu)；（2）每個流動層都視為均勻單孔隙度介質(zhì)。對于實際EGS熱儲，一般只有2～3分層[26]，本文將熱儲沿深度方向等分為10個水平流動層。

圖1 EGS地下部分和分層熱儲的物理模型Fig. 1 Physical model of EGS subsurface part and stratified reservoir

1.2 數(shù)學(xué)模型

對于EGS采熱過程的模擬，我們采用了前期自主開發(fā)的三維模型[13]。該模型為了模擬熱儲中的局部對流換熱，基于局部非熱平衡思想[30]，采用兩個能量方程來分別描述流體和巖石的溫度變化過程。模型的主要假設(shè)如下：

（1）單相流體流動；

（2）初始時刻，熱儲孔隙中充滿采熱流體；

（3）忽略循環(huán)流體的損失量（fluid loss）；

（4）不考慮循環(huán)流體與巖石的化學(xué)反應(yīng)；

（5）忽略巖石熱應(yīng)變及流體壓力導(dǎo)致的熱儲孔隙度的變化。

模型的主控方程有連續(xù)性方程、動量方程、流體的能量方程和巖石的能量方程。

連續(xù)性方程：

動量方程：

流體的能量方程：

巖石的能量方程：

其中，ε為熱儲孔隙度，ρ為密度，cp為比熱容，u為流體表觀速度，P為流體壓力，μ為流體動力粘性系數(shù)，K為熱儲滲透率；熱儲中巖石和流體的有效導(dǎo)熱系數(shù)使用修正因子為1.5的Bruggeman關(guān)系式來確定，即和h代表熱儲中固體巖石和流體間的對流換熱系數(shù)，a為裂隙的比表面積，可依據(jù)平板裂隙結(jié)構(gòu)或圓管直通道來估測h和a的數(shù)值[17]，研究中取ha=10 W/(m3·K)。

2 算例設(shè)置

2.1 模型參數(shù)

雙井EGS的幾何參數(shù)如圖2所示。人工熱儲的體積為500 m × 500 m × 500 m的正方體，計算區(qū)域大小為1 000 m × 1 000 m × 1 500 m，為減少人為邊界條件設(shè)置對結(jié)果帶來的影響，熱儲周圍包覆有足夠體積的巖石。注入井和生產(chǎn)井均為0.2 m × 0.2 m的方形通道。模型的求解借助Fluent通用計算流體力學(xué)軟件提供的求解平臺，通過其靈活的UDF（User Defined Functions）進行編程并數(shù)值求解。因幾何對稱，僅模擬一半的區(qū)域，數(shù)值網(wǎng)格劃分見圖 2，網(wǎng)格數(shù)約為250 000，該網(wǎng)格系統(tǒng)已經(jīng)進行過網(wǎng)格獨立性驗算，可以保證足夠的計算精度。

熱儲中心溫度為 180℃，當?shù)氐牡販靥荻葹?.04℃/m，熱儲大約位于地表下4 000 m位置，巖石邊界為定溫邊界。注入流體（水）的溫度為70℃。初始時，熱儲中裂隙流體與當?shù)貛r石溫度相同，注入井和生產(chǎn)井中都充滿70℃的流體。所有與流體接觸的壁面均為非滑移邊界，注入井入口定質(zhì)量流量，采出井出口定壓力。巖石和流體的熱物性參數(shù)設(shè)定為常數(shù)，如表1列示。

圖2 雙井EGS幾何參數(shù)和數(shù)值網(wǎng)格示意Fig. 2 Geometrical dimensions and numerical system of the considered doublet EGS

表1 流體和巖體的熱物性參數(shù)Table 1 Thermophysical properties of fluid and rock

2.2 熱儲構(gòu)建

考慮三個不同孔隙度分布的分層熱儲，包括兩個不均勻熱儲（HER1，HER2）和一個均勻熱儲（HOR）。三個熱儲的平均孔隙度一樣，都為0.01，熱儲孔隙度沿深度方向的分布如圖3所示。

熱儲不同層滲透率的確定依據(jù)Verma-Pruess模型[31]：

其中，εc表示臨界孔隙度，εc=0.8ε0；基準孔隙度ε0和滲透率K0分別取為0.01和1 × 10-14m2；增長指數(shù)n取為3。

圖3 孔隙度沿深度方向的變化Fig. 3 Porosity profile along the depth direction

3 分析與討論

為方便討論EGS的采熱性能，我們引入整體采熱率和局部采熱率的概念[16]。整體采熱率γ(t)（Heat extraction ratio）指在時刻t，EGS已經(jīng)開采的熱能與熱儲中可利用總能量的比值。局部采熱率γ1(t)（Local heat extraction ratio）指在時刻t，熱儲巖石局部位置已經(jīng)開采的熱量與當?shù)乜衫每偰芰康谋戎?。兩者定義式分別為：

其中，Ts和Ts(t)分別表示巖石的初始溫度和t時刻的溫度，Tf,in和Tf,out分別表示注入井和生產(chǎn)井流體的溫度，qv為流體工質(zhì)的循環(huán)體積流量。

文章中判斷熱穿透（Thermal breakthrough）出現(xiàn)的標準為開采過程中生產(chǎn)井周圍區(qū)域的局部采熱率達到0.01，定義EGS使用年限為生產(chǎn)井出口流體溫度下降10 K時EGS所經(jīng)歷的運行時間。

3.1 熱儲非均質(zhì)性的影響

3.1.1 采熱性能

圖4顯示了在定循環(huán)流量下（75 kg/s）3組不同熱儲EGS生產(chǎn)溫度和采熱率隨時間的變化情況，并示出了當EGS到達使用年限（出口井溫降為10 K）時的采熱率。

圖4 熱儲非均質(zhì)性對生產(chǎn)井溫度和采熱率的影響Fig. 4 The influence of reservoir heterogeneity on production temperature and heat extraction ratio

不同熱儲EGS對應(yīng)的曲線顯示出相同的變化趨勢，在采熱初期，出口井流體溫度會維持在某一較高采熱溫度，隨后會隨EGS的運行而逐漸降低。然而，不同熱儲EGS生產(chǎn)溫度變化速度明顯不同，相比于HOR EGS，HER1 EGS和HER2 EGS生產(chǎn)溫度開始下降的時間較早，并且溫度的下降也更加迅速，使后兩者更早地達到使用年限，導(dǎo)致所開采的熱總量較低。

3.1.2 滲流場和溫度場

為了進一步尋找不均勻熱儲采熱性能較差的原因，選取HOR EGS和HER1 EGS，對這兩者的滲流場和溫度場進行分析比較。

因為熱儲滲透率極低，井中的流阻相對熱儲內(nèi)Darcy阻力可以忽略[16]，流體在熱儲中沿深度方向的流速相對于水平方向（X，坐標系參見圖2）流速可以忽略。圖5示出了在75 kg/s工質(zhì)流量時，HOR和 HER1熱儲內(nèi)的 X方向流體滲流速度分布。在HER熱儲內(nèi)，流體工質(zhì)在達西阻力小的高滲透性區(qū)域形成明顯的優(yōu)勢流動。

圖5 定流量下（75 kg/s）HOR與HER1內(nèi)的滲流場Fig. 5 Seepage flow field in HOR and HER1 (75 kg/s flow rate)

圖6截取了HOR和HER1在5.357年和10.702年的熱量開采情況：a～d為熱儲溫度分布情況；e～h為熱儲局部采熱率分布情況，其中三個等值面的局部采熱率分別為0.99、0.20和0.01。

局部采熱率為0.01的等值面可以近似表示熱開采過程在熱儲中到達的前沿位置，當該等值面到達生產(chǎn)井時，認為該區(qū)域已經(jīng)出現(xiàn)熱穿透現(xiàn)象。如圖6所示，HER1中的優(yōu)勢流動區(qū)域熱開采速度明顯高于其他區(qū)域。在5.357年時，HER1中的優(yōu)勢流動區(qū)域已經(jīng)發(fā)生熱穿透（圖6f），導(dǎo)致HER1 EGS生產(chǎn)井溫度下降（參見圖4）。與之相比熱開采在HOR 更均勻進行，當開采時間進行到10.702年時，HOR才出現(xiàn)熱穿透現(xiàn)象（圖 6g）。另一方面，出現(xiàn)熱穿透之后，生產(chǎn)井附近的熱儲區(qū)域溫度會迅速降低，導(dǎo)致生產(chǎn)井溫度進一步下降，縮短了系統(tǒng)的使用年限。從圖6h中可以看出，當開采進行到10.702年時，在HER1的優(yōu)勢流動區(qū)域中，局部采熱率為0.99的等值面已經(jīng)迫近生產(chǎn)井，這時生產(chǎn)井出口流體溫降達到11.07℃，HER1 EGS已經(jīng)超過了使用年限。此時，HER1熱儲中部分區(qū)域的熱開采過程還未充分進行，導(dǎo)致了熱儲利用率較低。

圖6 HOR與HER1內(nèi)巖石溫度場（a～d）和當?shù)夭蔁崧史植记闆r（e～h）Fig. 6 Distribution of rock temperature (a～d) and local heat extraction ratio (e～h) in two reservoirs (HOR and HER1)

3.2 定量化評估

從以上分析中可以看出，熱儲在豎直方向非均質(zhì)性會引起系統(tǒng)采熱性能的下降。為了定量評估熱儲的非均質(zhì)性對系統(tǒng)采熱性能的影響，本文從以下三個方面展開工作：（1）分析熱儲結(jié)構(gòu)特性對采熱性能的影響，歸納總結(jié)適合描述不均勻熱儲結(jié)構(gòu)特征的統(tǒng)計參數(shù)；（2）比較HER1、HER2和HOR三種熱儲EGS在多種運行條件下的模擬結(jié)果，尋找能夠關(guān)聯(lián) EGS采熱性能和熱儲非均質(zhì)性的指標因子；（3）模擬更多具有不同熱儲結(jié)構(gòu)EGS的長期運行，根據(jù)模擬結(jié)果擬合出EGS采熱性能與熱儲非均質(zhì)性之間的函數(shù)關(guān)系。

3.2.1 熱儲不均勻度α

我們定義 EGS熱儲不均勻度α（Reservoir heterogeneity）以定量化描述熱儲的分層非均質(zhì)性，要求該參數(shù)能夠充分地體現(xiàn)熱儲非均質(zhì)性對EGS采熱過程的影響。從3.1.2的分析中可以看出，熱儲中流量分布的不均勻是引起EGS使用效率下降的關(guān)鍵因素。而熱儲中各個流動區(qū)域的流量的方差，是評價熱儲中流量分布不均勻程度最直觀的統(tǒng)計參數(shù)。因此我們采用類似方差的形式，用流量系數(shù)的m階偶數(shù)中心距定義熱儲不均勻度α來描述熱儲的非均質(zhì)性。其表達式為：

其中，N為分層熱儲模型的層數(shù)（本文中N=10），m為中心距的階次（m可取值2、4、6、8等）。流量系數(shù)σi為分層熱儲中第i個流動區(qū)域中的流量占總流量的比例（為所有流動層流量系數(shù)的平均值），其定義式為：

上式中，q為系統(tǒng)的總流量，qi為第i個流動區(qū)域中的流量。

在所考慮的分層熱儲模型中，流量系數(shù)可以進一步由每個流動區(qū)域的滲透率表示。流體工質(zhì)在井筒中因流動的粘性力產(chǎn)生的壓降要遠遠小于其在熱儲中由達西阻力導(dǎo)致的壓降，在熱儲滲透率小于10-12m-2時EGS的井筒中可近似看做定壓力[16]。所以圖1中的各個流動層可視為在兩井之間相互“并聯(lián)”，由達西定理[32]可知：

式中，Ai和A分別表示第i個流動區(qū)域在兩井間的等效截面積和熱儲區(qū)域的等效截面積（由于考慮的模型被均分為10個流動區(qū)域，因此Ai=A/10）；Ki和K分別表示第i個流動區(qū)域滲透率和熱儲區(qū)域的宏觀滲透率；μ為流體動力粘性系數(shù)；ΔP為入口井到生產(chǎn)井的壓降。

結(jié)合式（10）和式（11），第i個流動層的流量系數(shù)σi可直接由滲透率計算得出：

由式（8）和式（12）我們可以得到由熱儲各個流動層的滲透率計算得到熱儲的不均勻度α：

3.2.2 相對采熱率η

我們的目標是建立熱儲的非均質(zhì)性與EGS采熱性能的定量關(guān)系。為了增強該結(jié)果的通用性，我們希望所采用的評價EGS采熱性能的指標因子不僅能夠區(qū)分熱儲非均質(zhì)性，還要不受其他運行條件（流體流量、熱儲尺寸）的影響。為此我們比較了HER1、HER2和HOR三種熱儲EGS在多種運行條件下的模擬結(jié)果。

圖7顯示了三個熱儲（HOR、HER1、HER2）在不同流量下的采熱情況。圖7a顯示了EGS采熱率γ 隨生產(chǎn)井流體溫降δ（Production temperature decline）的變化，可以看出在相同流量下，均勻熱儲（HOR）的采熱率都明顯高于不均勻熱儲（HER1和HER2）。另一方面，EGS熱儲在熱開采過程中存在熱補償現(xiàn)象，熱儲會不斷地從周圍巖石中吸收熱量，循環(huán)流量較小時EGS熱開采過程進行緩慢，使熱儲能夠更充分地吸收周圍巖石的熱量，增加了EGS采熱總量，使系統(tǒng)采熱率發(fā)生改變：在熱儲HER1中，當流量從75 kg/s變?yōu)?7.5 kg/s、生產(chǎn)井流體溫降為10℃時的采熱率從0.295變?yōu)?.341，變化率為15.6%?？梢钥闯?，EGS采熱率雖然能夠區(qū)分熱儲的非均質(zhì)性，但循環(huán)流量的變化對EGS采熱率也有很大的影響。

為了突出熱儲非均質(zhì)性的影響，我們引入了相對采熱率η（Relative heat extraction ratio）。其定義為在相同的運行條件下生產(chǎn)井流體達到相同的溫降時系統(tǒng)的采熱率γ與均勻熱儲采熱率γ0的比值，表達式為：

如圖7b中所示，采熱曲線明顯地分成了三組。與圖7a相比，流量因素的影響幾乎完全消除。同樣，在熱儲HER1中，當流量從75 kg/s變?yōu)?7.5 kg/s、溫降為10 K時的相對采熱率從0.666變?yōu)?.686，變化率為3.0%。熱儲非均質(zhì)性成為影響系統(tǒng)相對采熱率的主導(dǎo)性因素。

圖7 EGS在不同流量下的采熱率（a）和相對采熱率（b）隨生產(chǎn)井溫降的變化Fig. 7 Heat extraction ration (a) and relative heat extraction ratio (b) versus production temperature decline for the EGS cases

對于不同的EGS，熱儲厚度、兩井間距離和熱儲區(qū)域?qū)挾瓤赡軙兴煌?。為了驗證相對采熱率對這些因素敏感性，模擬了當HER1熱儲模型在各個方向上的尺寸發(fā)生變化時的運行情況?？紤]到當熱儲長度縮短時，由于工質(zhì)在熱儲中流動距離的縮短，熱儲熱量的開采速度會相對變快。這與系統(tǒng)運行流量增大時的情況相似。為了避免計算結(jié)果與之重復(fù)，我們在設(shè)置熱儲長度縮小的算例時，以相同的比例減小了系統(tǒng)的運行流量。結(jié)果如圖8所示。圖中l(wèi)x、ly、lz分別表示立方體熱儲在x、y、z三個方向上的邊長，當立方體邊長改變時熱儲內(nèi)部尺寸也以相同比例變化。

可以看出，當熱儲尺寸發(fā)生變化時不均勻熱儲的相對采熱率并無明顯變化。從局部放大圖中看到，熱儲長度lx的減小幾乎不引起相對采熱率的改變，說明兩井距離并不會影響系統(tǒng)的相對采熱率。熱儲高度ly的減小會導(dǎo)致相對采熱率的小幅上漲，而熱儲寬度lz的減少會導(dǎo)致相對采熱率的小幅下降，但變化幅度十分微弱。說明流動層厚度和熱儲寬度的改變都不會對EGS相對采熱率產(chǎn)生明顯影響。

圖8 不同尺寸HER1 EGS相對采熱率隨生產(chǎn)井溫降的變化Fig. 8 Relative heat extraction ratio versus production temperature decline for HER1 EGSs of differing sizes

通過以上分析，我們看到相對采熱率η不僅可以有效表征熱儲非均質(zhì)性對EGS采熱的影響，而且η-α關(guān)系受流體循環(huán)流量和熱儲尺寸的影響很小。采用相對采熱率η作為描述 EGS采熱性能的指示因子，可以使研究結(jié)果更具通用性。

3.2.3 η與α擬合關(guān)系

為探求普遍情形下EGS中相對采熱率η和熱儲不均勻度α之間的關(guān)系，形成系統(tǒng)的、定量的研究結(jié)論。本文通過分析多個不均勻熱儲模擬結(jié)果，嘗試總結(jié)出在一定情形下具有普適意義的關(guān)系式。為此，我們構(gòu)建了30個不同的熱儲，每個不均勻熱儲都劃分為10個孔隙度不同的熱儲層（參考HER1和HER2），而每個熱儲層的孔隙度都利用Matlab中高斯分布隨機數(shù)生成器（Normrnd函數(shù)）隨機賦值，并保證10個熱儲層的平均孔隙度為0.01。此外，為了對不均勻熱儲的不均勻程度加以區(qū)分，30個熱儲孔隙標準差（Standard Deviation, STD）在0.01～1 × 10-3范圍內(nèi)變化。模擬時循環(huán)流量統(tǒng)一設(shè)為75 kg/s，熱儲都是邊長為500 m的正方體，并采用系統(tǒng)達到使用年限（生產(chǎn)井溫降δ=10 K）時的相對采熱率作為評價指標。

通過對模擬結(jié)果進行擬合處理，發(fā)現(xiàn)m=4時的熱儲不均勻度α與系統(tǒng)相對采熱率η之間存在很好的關(guān)聯(lián)性。如圖9所示，當用熱儲孔隙度的標準差作為橫坐標時，標準差大于0.3 × 10-3的數(shù)據(jù)點分布分散，熱儲孔隙度的標準差與系統(tǒng)的相對采熱率并沒有很大的關(guān)聯(lián)性；當使用熱儲不均勻度α（m=4）作為橫坐標時，數(shù)據(jù)點的離散程度明顯降低，數(shù)據(jù)點可以使用一個雙曲正切函數(shù)進行擬合：

該擬合結(jié)果的平均相對誤差為 2.24%，最大相對誤差為 7.45%（該相對誤差為擬合誤差與模擬結(jié)果值的比值）。

圖9 30個算例的模擬結(jié)果與擬合曲線（考慮生產(chǎn)井溫降為10 K時的相對采熱率，取m=4計算熱儲不均勻）Fig. 9 η-αrelationship based on data from 30 simulations (ηvalues taken at 10K production temperature decline, andαvalues calculated withm=4)

從擬合結(jié)果可以看出，當熱儲不均勻度α的量級為 10-6時，系統(tǒng)的相對采熱率η低于 0.9；當α進一步增大時，η出現(xiàn)顯著下降，此時熱儲的非均質(zhì)性對系統(tǒng)的采熱性能造成了較大影響。

4 總 結(jié)

考慮到熱儲在豎直方向的非均質(zhì)性，本文將不均勻熱儲假設(shè)為分層等效多孔介質(zhì)。研究發(fā)現(xiàn)，分層不均勻熱儲中各個流動層流體流量的差別是導(dǎo)致系統(tǒng)采熱性能下降的最主要原因。

在分層熱儲模型基礎(chǔ)上，本文引入一種全新的熱儲非均質(zhì)性描述和評價方法：由熱儲的滲透率分布情況計算出熱儲不均勻度α，并使用相對采熱率η來表示不均勻熱儲采熱性能的下降。

本文通過隨機賦值方法構(gòu)建了多個不同的熱儲結(jié)構(gòu)，并模擬了它們各自的長期運行情況?；谶@些模擬結(jié)果，發(fā)現(xiàn)使用各個流動層流量系數(shù)的四階中心距α能夠很好地關(guān)聯(lián)熱儲非均質(zhì)性與EGS采熱性能，并給出了η和α的擬合關(guān)系式。

根據(jù)該關(guān)系式，可以直接利用不均勻熱儲滲透率的分布情況對具有分層熱儲結(jié)構(gòu)的EGS采熱性能進行定量預(yù)測。不均勻熱儲在不同的運行流量、熱儲尺寸下的運行情況的模擬結(jié)果顯示該η與α關(guān)系式受流體流量、熱儲尺寸影響較小，這表明該關(guān)系式具有良好的通用性。

[1] Tester J W, Anderson B J, Batchelor A S, et al. The future of geothermal energy[R]. Massachusetts Institute of Technology, 2006.

[2] 藺文靜, 劉志明, 馬峰, 等. 我國陸區(qū)干熱巖資源潛力估算[J]. 地球?qū)W報, 2012, 33(5): 807-811.

[3] 汪集暘, 胡圣標, 龐忠和, 等. 中國大陸干熱巖地?zé)豳Y源潛力評估[J]. 地?zé)崮? 2013, (3): 3-10.

[4] 許天福, 張延軍, 曾昭發(fā), 等. 增強型地?zé)嵯到y(tǒng)(干熱巖)開發(fā)技術(shù)進展[J]. 科技導(dǎo)報, 2012, 30(32): 42-45.

[5] 郭劍, 陳繼良, 曹文炅, 等. 增強型地?zé)嵯到y(tǒng)研究綜述[J]. 電力建設(shè), 35(4): 10-24.

[6] Breede K, Dzebisashvili K, Liu X, et al. A systematic review of enhanced (or engineered) geothermal systems: past, present and future[J]. Geothermal Energy, 2013, 1(1): 4.

[7] McClure M W. Fracture stimulation in enhanced geothermal systems[D]. California: Stanford University, 2009.

[8] Dezayes C, Genter A, Hooijkaas G R. Deep-seated geology and fracture system of the EGS Soultz reservoir (France) based on recent 5km depth boreholes[C]// Proceedings World Geothermal Congress 2005, Antalya, Turkey, 24-29 April, 2005.

[9] Genter A, Evans K, Cuenot N, et al. Contribution of the exploration of deep crystalline fractured reservoir of Soultz to the knowledge of enhanced geothermal systems (EGS)[J]. Comptes Rendus Geoscience, 2010, 342(7): 502-516.

[10] Jung A, Fenwick D H, Caers J. Training image-based scenario modeling of fractured reservoirs for flow uncertainty quantification[J]. Computational Geosciences, 2013, 17(6): 1015-1031.

[11] Berkowitz B. Characterizing flow and transport in fractured geological media: A review[J]. Advances in water resources, 2002, 25(8): 861-884.

[12] Dershowitz W S, La Pointe P R, Doe T W. Advances in discrete fracture network modeling[C]//Proceedings of the US EPA/NGWA Fractured Rock Conference, Portland. 2004: 882-894.

[13] Jiang F, Luo L, Chen J. A novel three-dimensional transient model for subsurface heat exchange in enhanced geothermal systems[J]. International Communications in Heat and Mass Transfer, 2013, 41: 57-62.

[14] 陳繼良, 蔣方明. 增強型地?zé)嵯到y(tǒng)熱開采性能的數(shù)值模擬分析[J]. 可再生能源, 2013, 31(012): 111-117.

[15] 陳繼良, 羅良, 蔣方明. 熱儲周圍巖石熱補償對增強型地?zé)嵯到y(tǒng)采熱過程的影響[J]. 計算物理, 2013, 30(6): 862-870.

[16] 陳繼良, 蔣方明, 羅良. 增強型地?zé)嵯到y(tǒng)地下滲流場的模擬分析[J]. 計算物理, 2013, 30(6): 871-878.

[17] Arbogast T, Douglas, Jr J, Hornung U. Derivation of the double porosity model of single phase flow via homogenization theory[J]. SIAM Journal on Mathematical Analysis, 1990, 21(4): 823-836.

[18] Gelet R, Loret B, Khalili N. A thermo-hydro-mechanical coupled model in local thermal non-equilibrium for fractured HDR reservoir with double porosity[J]. Journal of Geophysical Research: Solid Earth (1978-2012), 2012, 117(B7): B07205-B07228.

[19] Pruess K. The TOUGH codes-A family of simulation tools for multiphase flow and transport processes in permeable media[J]. Vadose Zone Journal, 2004, 3(3): 738-746.

[20] Xu T, Spycher N, Sonnenthal E, et al. TOUGHREACT Version 2.0: A simulator for subsurface reactive transport under non-isothermal multiphase flow conditions[J]. Computers & Geosciences, 2011, 37(6): 763-774.

[21] Kalinina E, McKenna S A, Hadgu T, et al. Analysis of the effects of heterogeneity on heat extraction in an EGS represented with the continuum fracture model[C]// Proceedings, 37th Workshop on Geothermal Res. Eng., Stanford Univ., Stanford, CA, Jan. 2012.

[22] Shaik A R, Rahman S S, Tran N H, et al. Numerical simulation of Fluid-Rock coupling heat transfer in naturally fractured geothermal system[J]. Applied Thermal Engineering, 2011, 31(10): 1600-1606.

[23] Vogt C, Marquart G, Kosack C, et al. Estimating the permeability distribution and its uncertainty at the EGS demonstration reservoir Soultz-sous-Forêts using the ensemble Kalman filter[J]. Water Resources Research, 2012, 48(8): W08517-W08531.

[24] Vogt C, Kosack C, Marquart G. Stochastic inversion of the tracer experiment of the enhanced geothermal system demonstration reservoir in Soultz-sous-Forêts-Revealing pathways and estimating permeability distribution[J]. Geothermics, 2012, 42: 1-12.

[25] Fourar M. Characterization of heterogeneities at the core-scale using the equivalent stratified porous medium approach[C]//SCA International Symposium. Trondheim, Norway. 2006.

[26] Radilla G, Sausse J, Sanjuan B, et al. Interpreting tracer tests in the enhanced geothermal system (EGS) of Soultz-sous-Forêts using the equivalent stratified medium approach[J]. Geothermics, 2012, 44: 43-51.

[27] Bl?cher M G, Zimmermann G, Moeck I, et al. 3D numerical modeling of hydrothermal processes during the lifetime of a deep geothermal reservoir[J]. Geofluids, 2010, 10(3): 406-421.

[28] Luo F, Xu R N, Jiang P X. Numerical investigation of fluid flow and heat transfer in a doublet enhanced geothermal system with CO2as the working fluid (CO2-EGS)[J]. Energy, 2014, 64: 307-322.

[29] Luo F, Xu R N, Jiang P X. Numerical investigation of the influence of vertical permeability heterogeneity in stratified formation and of injection/production well perforation placement on CO2geological storage with enhanced CH4recovery[J]. Applied Energy, 2013, 102: 1314-1323.

[30] He L W, Jin Z H. A local thermal nonequilibrium poroelastic theory for fluid saturated porous media[J]. Journal of Thermal Stresses, 2010, 33(8): 799-813.

[31] Lai K H, Chen J S, Liu C W, et al. Effect of permeability–porosity functions on simulated morphological evolution of a chemical dissolution front[J]. Hydrological Processes, 2014, 28(1): 16-24.

[32] 孔祥言. 高等滲流力學(xué)[M]. 合肥: 中國科學(xué)技術(shù)大學(xué)出版社, 1999. 30-37.

Numerical Study of EGS Heat Extraction from Stratified Heat Reservoirs

HUANG Wen-bo, CHEN Ji-liang, JIANG Fang-ming
(Laboratory of Advanced Energy System, CAS Key Laboratory of Renewable Energy, Guangzhou Institute of Energy Conversion, Chinese Academy of Sciences, Guangzhou 510640, China)

Strong heterogeneity exists in enhanced geothermal systems (EGS) artificial heat reservoir due to factors such as different rock formations, the pre-existence of natural fractures and the uncertainty of stimulation process. Exploring the influence of reservoir heterogeneity on EGS heat extraction is of great importance for predicting and comprehensively evaluating the performance of EGS. In this paper, we take the reservoir as an equivalent stratified porous medium and employ the previously self-developed numerical model to conduct a series of simulations of the long-term heat extraction process of EGSs with different reservoir porosity distributions. Results indicate that the heat extraction performance is significantly affected by the reservoir heterogeneity. To facilitate the analysis and evaluation, we derive a formula to quantitatively describe the heterogeneity of the reservoir based on a detailed analysis to the effects of the fluid seepage flow field on the heat extraction performance. Further, we obtain a quantitative relationship of EGS heat extraction performance versus reservoir heterogeneity by fitting the calculated results from 30 cases, in which the stratified reservoirs are randomly assigned with differing porosity distributions.

enhanced geothermal systems (EGS); numerical simulation; equivalent stratified porous medium; heterogeneity

TK529；P314

A

10.3969/j.issn.2095-560X.2014.04.009

2095-560X（2014）04-0295-10

黃文博（1990-），男，碩士研究生，主要從事增強型地?zé)嵯到y(tǒng)地下物理過程的數(shù)值模擬研究。

2014-05-23

2014-07-09

國家高技術(shù)發(fā)展863計劃項目（2012AA052802）；中國科學(xué)院“百人計劃”項目

? 通信作者：蔣方明，E-mail：jiangfm@ms.giec.ac.cn

蔣方明（1973-），男，博士，研究員，博士生導(dǎo)師，中國科學(xué)院廣州能源研究所先進能源系統(tǒng)實驗室主任。2001年博士畢業(yè)后曾先后工作于德國的 IMM 公司、葡萄牙的阿維諾（Aveiro）大學(xué)、美國的賓州州立大學(xué)，2011年4月獲中國科學(xué)院“百人計劃”引進海外杰出人才擇優(yōu)支持。目前主要從事電化學(xué)能量/動力系統(tǒng)、增強型地?zé)嵯到y(tǒng)、微熱流體系統(tǒng)、以及高效節(jié)能技術(shù)/產(chǎn)品等研發(fā)工作。