郑州市主城区与东部新城区地热资源数值模拟研究

王 宽,周明泽,张一博,马孟科

(1.河南省地质局地质灾害防治中心,河南 郑州 450045;2.河南省自然资源科技创新中心(豫北地热能清洁能源研究),河南 郑州 450045)

随着国家经济的不断增长,能源需求方面也不断加大,人均能源消耗量在增加,但是要求碳排放量要大幅降低。这就要求清洁能源发挥强大作用,亟须能源结构转型。地热能资源是清洁环保的新型可再生能源,在我国未来能源结构调整中发挥重要作用[1-3]。随着我国城镇化、工业化进程加快和经济社会的迅速发展,资源与环境矛盾日趋突出,地热清洁能源的开发利用越来越得到国家和社会的重视,依托深部地热供暖技术扩大地热能资源开发利用规模也成为节能减排、优化能源结构、减少空气污染、提高人民生活品位的重要措施之一[4-7]。

郑州市以往开展的地热资源勘查评价工作虽然取得了比较丰富的成果,但大多数仅停留在地热地质研究及资源评价层面,在地热资源勘查手段、开采工艺、开发利用研究方面较少。总体上郑州市地热开发利用缺乏科学的规划指导,导致部分地区开采分布集中,而开采量也没有经过科学规划,开采后不进行地下水的补充,最终造成地下热水资源水位和水温下降,影响开采后的收益,同时也破坏了资源存在的环境,影响了资源的可持续利用。

本文选取郑州市主城区及东部新城区为研究对象,建立地下水热耦合数值模拟模型,将地热资源开采对地下水环境产生的影响进行深入研究,并为合理的开采回灌井提供设计依据。

1 地热地质条件

1.1 地热系统分区

本次研究区地热类型主要划分为2类,呈层状分布的盆地型地热及主要受断裂构造控制呈带状分布地热。基于研究区构造特征,将研究区划分为3个地热区:郑州断阶地热区、中牟凹陷地热区和嵩箕台隆地热区(图1)。

图1 地热系统分区Fig.1 Geothermal system zoning map

1.2 地温场特征

(1)地温梯度分布规律。根据研究区地热井调查资料,地温梯度变化范围2.54～4.87 ℃/hm,大部分地区地温梯度2.5～3.5 ℃/hm,局部存在地温异常。主城区平均地温梯度3.07 ℃/hm,东部新城区平均地温梯度3.27 ℃/hm。

(2)热储温度分布规律。以计算所得新生界地温梯度为基础,结合不同地层埋深分别计算了埋深1 500、2 000 m处的地温和新近系明化镇和馆陶组底板埋深的地温。1 500 m埋深地温变化范围50～78 ℃。其中,郑州断阶地热区温度在50～60 ℃,开封凹陷地热区温度在60～68 ℃,嵩箕台隆地热区温度在55～78 ℃。2 000 m埋深地温变化范围60～94 ℃。其中,郑州断阶地热区温度在60～75 ℃,开封凹陷地热区温度在75～85 ℃,嵩箕台隆地热区温度在66～94 ℃。新近系明化镇底板埋深的地温整体从西向东北逐渐增加,变化范围30～68 ℃。其中,郑州断阶地热区温度在30～50 ℃,开封凹陷地热区温度在50～68 ℃。新近系明化镇和馆陶组底板埋深的地温整体从西向东北逐渐增加,变化范围50～100 ℃。其中,郑州断阶地热区的温度在50～70 ℃,开封凹陷地热区温度在70～100 ℃。

1.3 热储特征及埋藏条件

(1)新近系明化镇组热储层。新近系明化镇组(N2m)热储层分布于工作区尖岗断层以北大部地区,主要分布在郑州断阶地热区和中牟凹陷地热区内。顶板埋深160～380 m。底板埋深自西南向东北逐渐加深,郑州断阶底板埋深一般在400～1 200 m,中牟凹陷底板埋深一般在1 100～1 800 m。含水层以中细砂为主,共有10余层,下部微胶结,属于半固结碎屑岩类孔隙裂隙含水层组,总厚度56～744 m,平均厚308 m。据区内钻孔统计,砂层厚度比平均46%左右。该储层大部分地区地热井地热流体单位产量大于50 m3/(d·m),为产流能力强区;郑州市西北部石佛一带—河南省经贸委一带、郑州工程机械厂—省送变电公司—郑州市南郊五里堡村—香水皇宫假日酒店一带产流能力中等,单位产量5～50 m3/(d·m);河南省华垦实业公司一带地热井地热流体单位产量弱,小于5 m3/(d·m)。该储层水温一般为25～40 ℃,为温水储层。

(2)新近系馆陶组热储层。新近系馆陶组(N1g)热储层分布受老鸦陈断层控制,主要分布在郑州断阶和中牟凹陷内。其顶板埋深大于800 m,底板埋深自西南郑州断阶向东北中牟凹陷逐渐加深,郑州断阶底板埋深一般在800～2 000 m,中牟凹陷底板埋深一般在1 300～2 800 m。热储层岩性下段为棕红、棕色半胶结泥岩与细中砂、中砂互层;中段为棕红、棕黄色半胶结泥岩与灰黄、灰白、黄白色中砂、中粗砂互层;上段为灰绿、灰棕、棕红色半胶结泥岩夹黄白色细砂、中砂,厚200～1 000 m。储水介质为半胶结的细砂、中细砂、砂砾石,下部微胶结,属于半固结碎屑岩类孔隙裂隙含水层组,分为8～10层,总厚度50～173 m,平均厚112 m。据区内钻孔统计,砂层厚度比平均42%左右。地热产流能力大部分地区为中等,地热井地热流体单位产量5～50 m3/(d·m)。在省老干部活动中心—郑东新区馨悦苑小区一带产流能落弱,地热井地热流体单位产量小于5 m3/(d·m);市区北部柳林村委会一带产流能力强,地热井地热流体单位产量大于50 m3/(d·m)。

2 地下水热耦合数值模拟

运用TOUGH2—EOS1模块建立郑州市主城区和东部新城区地下水热耦合数值模拟模型,根据现状条件下的地下水热信息资料对模型进行拟合验证,确保模型的可靠性[8-10],最后对模型重要参数进行了敏感性分析。

2.1 数值模型

(1)初始压力场、温度场。模型Ⅰ:模型设为稳定流,初始压力场—由于TOUGH模型计算采用地下水流体压力数据,则在计算过程中,根据初始时刻对应的钻孔中的水位观测资料,将地下水位换算位热储层压力值为模型赋值压力场。系统的初始压力分布通常由重力平衡确定:①通过插值得出重力平衡的初始计算条件;②对定水头网格块设定fixed state作为定压条件;③模拟长时间直至稳定,得到重力平衡下初始压力分布。模型初始流场如图2所示。

图2 模型初始流场Fig.2 Initial flow field of model

初始温度场—根据研究区地温场特征,确定恒温带深度27 m,恒温带温度17 ℃,地温梯度根据地热能赋存条件,参考研究区地温梯度取3 ℃/hm,为增温带温度场进行赋值,地层温度计算如式(1)。

T1=T0+(S1-S0)×G/100

(1)

式中,T1为计算深度地温;T0为恒温带温度;S1为计算深度;S0为恒温带深度;G为地温梯度。

模型Ⅱ:上述初始模型运行稳定,将运行结果作为下阶段模型的初始条件。模型初始流场如图3所示。

图3 模型初始温度场Fig.2 Initial heat field of model

(2)边界条件。水平方向上,东、西边界处理为流量边界,其中西边界为流入边界,东边界为流出边界,利用郑州市新近系热储层等水位线图,计算出模拟区的水力梯度大约为1/1 000,流入流出量通过达西定律计算得到;南北边界处理为隔水边界。

2.2 参数选取及分区

(1)水文地质参数及热物性参数。水文地质参数及热储物理性质参数包括:岩石密度、孔隙度、比热和热导率等。

本次工作根据各热储的埋藏分布和控热构造分布特征,结合收集的试验数据分析整理得到。

其中岩石密度、孔隙度参考钻孔岩性取样检测报告,渗透率参考工作区钻孔渗透试验数据分析得到,研究区参数分区设置及参数值见表1、表2。

表1 研究区参数分层设置及参数值(热储Nm)Tab.1 Hierarchical setting and parameter value of research area(thermal storage Nm)

表2 研究区参数分层设置及参数值(热储Ng)Tab.2 Hierarchical setting and parameter value of research area(thermal storage Ng)

(2)大地热流值。是地球内热在地表唯一可以量测的物理量,比其他地热参数更能确切地反映某个地区地温场的特点。其计算公式为:

q=-100Kr(dT/dz)

(2)

式中,q为大地热流,通常缩写为HFU(Heat Flow Unit);Kr为岩石导热率;dT/dz为地温梯度。

2.3 模型识别验证

2.3.1 地下水位识别验证

本次模型选取收集到研究区部分观测井平均水位和长期水位监测钻孔G4、G6、G7及地下水环境监测站CR306的水位数据进行识别验证。经模拟调参运行,拟合结果显示见表3,观测井的实测平均水位与模拟水位误差符合地热资源评价规范中要求,拟合结果较为理想,模拟结果可以反映地下水流场特征,可用于解释地下水流动模式。

表3 模拟水位与实测水位对比Tab.3 Comparison between simulated water level and measured water level

2.3.2 地下水温度场识别验证

根据研究区地热井温度监测数据,本次模拟假设模拟期内各地热井温度保持不变。钻孔温度拟合曲线如图4所示,误差符合地热资源评价相关规范的要求。研究区深层地温场温度随热储层埋深增加而相应升高,模拟地下水温度场符合实际情况。

图4 钻孔温度实测、模拟结果对比曲线Fig.4 Comparison curve between measured and simulated borehole temperature

2.4 均衡分析

模型现状开采条件下水均衡见表4。模型总补给量为38.1万m3,总排泄项量为180.7万m3,均衡差为-142.6万m3。郑州市现状开采条件下地下热水主要排泄项为人工开采,历来多为饮用、洗浴等生活用水以及供暖期取暖用水,供暖井大部分虽配有回灌井,但受回灌量与回灌技术限制,地热尾水回灌比例较小,研究区热储层地下水系统仍为负均衡状态,地下水位有持续降低的风险。

表4 水均衡Tab.4 Water balance

3 采灌模拟分析

综合郑州市地热资源发展规划及研究区地质勘察相关资料,分别对郑州市馆陶组及明化镇组热储层进行采灌方案设计并模拟分析不同方案对地下水位及温度的影响规律,最后对模型重要参数进行了敏感性分析[11-14]。

3.1 馆陶组模拟

3.1.1 馆陶组模拟工况方案设计

由于本次模型较大,为了更具体、更清晰分析地下热水抽回灌模拟的水位及温度变化规律,以建立的模型为基础,选取东部新城区馆陶组热储层作为地热井的采灌位置来设计不同模拟工况方案,小模型边长为5 km×5 km,底板标高为-2 200 m,馆陶组热储层开采位置为-1 300～ -2 200 m。以馆陶组热储层在“一抽一灌”及“一抽两灌”的不同采灌方式下,地热井开采及回灌只发生在一个供暖周期内(供暖120 d),剩余时间用于恢复水位及含水层热量,地热井采灌系统设计工作年限为30年,开采量及回灌量为60 m3/h,回灌井的回灌温度为15 ℃情况下,共设计模拟6种工况,见表5,分析采灌系统运行30年不发生热突破的合理井间距及回灌方式,并分析不同工况条件对采灌井水位及温度场的影响。

表5 模拟工况汇总Tab.5 Summary of simulated working conditions

3.1.2 馆陶组模拟结果分析

(1)热均衡分析。模型整体源汇项包含侧向水流的补给与排泄以及开采井的开采与回灌井的回灌,由于模拟工况均为100%回灌模拟,而侧向的补给与排泄量大致相等,故模型地下水流系统相对均衡。通过对模型热量结果文件整理分析,总体补给热量为2.03×1015J,侧向排泄热量仅为3.09×1012J,大地热流为主要热储补给来源,根据相关文献研究,当热储层深度越深,大地传导热流的补给越显著,地热井的温度升高也就越明显,所以温度随地热井深度的增加而逐渐升高。开采井开采的热量与模型侧向排泄及顶部向周围岩石传递热量为主要热量输出渠道。

(2)“开采井、回灌井”“热突破”模拟分析。馆陶组热储层各工况的模拟情况如图5所示(红色为回灌井,蓝色为开采井),各模拟工况供暖期末开采井温度变化情况如图6所示。在地热井运行期受回灌井回灌温度的影响,开采井温度会出现一定程度降低,发生热突破现象,开采井温度变化幅度受采灌井系统运行时间、不同的回灌温度以及研究区水文地质条件影响,结合本次研究内容以开采井的温度下降 0.3 ℃作为形成热突破的标准。由于热对流扩散现象,等温线主要分布在回灌井附近,沿地下水流向开采井方向扩展。“一抽一灌”工况下,采灌井间距150 m时在24年左右发生热突破,运行30年后开采井温度下降约0.5 ℃,采灌井间距200、250 m时均未发生热突破;“一抽两灌”工况下,采灌井间距150 m时运行30年后开采井温度下降约0.3 ℃,采灌井间距200、250 m时均未发生热突破。

图5 各工况条件下温度影响范围Fig.5 Temperature influence range under various working conditions

图6 各工况供暖期末开采井温度变化Fig.6 Temperature change of mining well at the end of heating period under various working conditions

3.1.3 采灌方案结果影响分析

根据模拟分析结果,进一步对井间距为150、200 m的模拟工况进行分析,对比其“一抽一灌”及“一抽两灌”条件下采灌区的渗流场及温度场变化情况,通过对比水位及温度的差异来探讨不同回灌方式及井间距对渗流场及温度场的影响。

150、200 m井间距“一抽一灌”、“一抽两灌”条件下供暖期末采灌区地下水位如图7所示。由图7可知,在采灌方式一致时,随着井间距增大,热储层供暖期末回灌井处水位呈现上升的趋势,开采井中心处的水位降深及水位漏斗面积增大。在井间距一致时,“一抽两灌”相较于“一抽一灌”条件下,热储层供暖期末回灌井与开采井中心处的水位及漏斗面积影响呈现减小的趋势。

图7 馆陶组供暖期末采灌区地下水位Fig.7 Groundwater level in the irrigation area at the end of heating period of Guantao Group

综上分析可知,本次模拟方案的地热井开采及回灌只发生在一个供暖周期内(供暖120 d),地热井采灌系统设计工作年限为30年,地热井开采量及回灌量为60 m3/h,回灌井回灌温度为15 ℃,在防止开采井发生热突破的前提下,馆陶组热储层“一抽一灌”工况方案采灌井的合理井间距的下限为200 m,“一抽两灌”工况方案采灌井合理井间距的下限为150 m。

3.2 明化镇组模拟

3.2.1 明化镇组模拟工况方案设计

选取主城区明化镇组热储层作为地热井的采灌位置来设计不同模拟工况方案,在模型回灌试验模拟中对采灌井位置处进行了网格加密,小模型边长为4 km×4 km,底板标高-620 m,馆陶组热储层开采位置为-200～ -600 m。以馆陶组热储层在“一抽一灌”及“一抽两灌”的不同采灌方式下,地热井开采及回灌仅在供暖期内进行(供暖时间为120 d),剩余时间用于恢复地热井的水位及含水层热量,地热井采灌系统设计工作年限为30年,开采量及回灌量均为60 m3/h,回灌井的回灌温度为15 ℃,设计不同模拟方案,分析采灌系统运行30年不发生热突破的合理井间距及回灌方式,并分析不同工况条件对采灌井水位及温度场的影响。

3.2.2 明化镇组模拟结果分析

(1)采灌井“热突破”模拟分析。本次研究明化镇组热储层模拟以开采井的温度下降 0.3 ℃作为形成热突破的标准。“一抽一灌”工况下,采灌井间距150 m时在22年左右发生热突破,运行30年后开采井温度下降约0.9 ℃,采灌井间距200 m运行30年后开采井温度下降约0.16 ℃,采灌井间距250 m运行30年后开采井温度下降约0.03 ℃,均未发生热突破;“一抽两灌”工况下,采灌井间距150 m时运行30年后开采井温度下降约0.35 ℃,采灌井间距200 m运行30年后开采井温度下降约0.05 ℃,采灌井间距250 m运行30年后开采井温度下降约0.02 ℃,均未发生热突破。

(2)采灌方案结果影响分析。将不同模拟工况进行研究分析,并对比“一抽一灌”及“一抽两灌”条件下采灌区的流场及温度场变化情况,通过讨论分析水位及温度的差异进而研究不同回灌方式及井间距对流场及温度场的影响。供暖期末采灌区地下水位如图8所示。

图8 明化镇组供暖期末采灌区地下水位Fig.8 Groundwater level in the irrigation area at the end of heating period of Minghua Town Group

从图8可以看出,明化镇组采灌模拟与馆陶组采灌模拟水位变化规律一致,随着井间距增大,热储层供暖期末回灌井处水位呈现上升的趋势,开采井中心处的水位降深及漏斗面积增大。在井间距一致时,“一抽两灌”相较于“一抽一灌”条件下,热储层供暖期末回灌井水位与开采井中心处的水位降落漏斗面积影响呈现减小的趋势,开采井发生热突破风险也相应减小。由于明化镇组热储层渗透率略高于馆陶组热储层,因此明化镇组水位变化大于馆陶组。

综上分析可知,明化镇组热储层在防止开采井发生热突破的前提下,“一抽一灌”工况方案采灌井的合理井间距的下限为250 m,“一抽两灌”工况方案采灌井合理井间距的下限为200 m,在实际工程中,选用“一抽两灌”开采方式可有效减少开采井中水位降深,减小降落漏斗面积,并且对地下水位及地下温度场影响更小。

3.3 参数敏感性分析

3.3.1 敏感性分析方法

为了研究各参数对含水层温度场、流场以及最佳采灌方案的影响程度,本节在前文地下水热耦合模型的基础上,以井间距200 m、“一抽一灌”的采灌方案为例,选取渗透率、岩石热传导系数以及不同回灌温度作为研究参数,以采灌井水位与开采井温度作为本次敏感性分析的模型输出结果进行参数敏感性分析。敏感度指数由式(3)确定:

(3)

式中,L为敏感度指数;ak为目标参数;m为分析组数;yi为分析结果。

3.3.2 敏感性分析结果

(1)渗透率(k)及岩石热传导系数(λ)。各参数增大50%及减小50%情况下的温度场影响范围如图9所示。由图9可知,渗透率及岩石热传导系数的变化引起的采灌井整体温度影响范围变化不大;其中,岩石热传导系数减小50%,模型运行30年后,回灌井处的圆形低温区面积有所增加,温度影响范围增大,开采井温度较参数不变时下降约0.035 ℃,岩石热传导系数增加50%时,模型运行30年后,回灌井处圆形低温区面积也相应减小,开采井温度较参数不变时上升约0.02 ℃;而渗透率的变化对采灌井的温度无明显影响。

图9 参数不同变幅下回灌温度影响范围Fig.9 Influence range of reinjection temperature under different parameter amplitudes

各参数不同变幅情况下采灌井水位统计见表6。由表6可知,岩石热传导系数对采灌区的水位影响较小,而渗透率改变造成采灌井水位有一定变幅。

表6 各参数不同变幅情况下采灌井中心水位Tab.6 Central water level of production and irrigation well under different amplitude of parameters

渗透率减小50%时,供暖期末回灌井中心水位较参数不变时上升1 m,水位上升区面积也有所增大,开采井中心水位较参数不变时降低0.5 m,且水位降落漏斗面积也相应增大;渗透率增加50%时,供暖期末回灌井中心水位较参数不变时降低0.34 m,水位上升面积也有所减小,开采井中心水位较参数不变时水位上升0.15 m,水位降落漏斗面积相应减小。

将不同渗透率K条件下采灌井水位变化展示如图10所示,由此分析可知随着渗透率的增大,开采井水位下降趋势逐渐增加,回灌井水位上升面积逐渐减小。

图10 供暖期末不同渗透率下采灌井间剖面线水位对比Fig.10 Comparison of profile water level between production and irrigation wells under different permeability at the end of heating period

(2)回灌温度。在实际工程中,地热水回灌温度根据地热水利用情况以及不同季节温度会有所不同,确定不同的回灌温度对采灌系统影响尤为重要,常见地热尾水回灌温度在10～50 ℃,冬季回灌温度低,夏季回灌温度高,因此模拟选取井间距200 m“一抽一灌”模拟工况,将模型的回灌温度由初始的15 ℃分别调整设置为30、45 ℃进行模拟,其余参数保持不变,分析得出不同回灌温度对采灌井的水位、温度场的影响范围。不同回灌温度下回灌井温度影响范围如图11所示。

图11 不同回灌温度下回灌井温度影响范围Fig.11 Influence range of reinjection well temperature under different reinjection temperatures

在不同回灌温度下,分别在回灌井处形成以15、30、45 ℃为中心的圆形低温回灌区,随着回灌温度的增加,中部低温区影响范围减小,200 m井间距范围、3种回灌温度条件下均未发生热突破,回灌温度15 ℃与30 ℃时开采井温度几乎无变化,回灌温度45 ℃时,开采井温度略微升高。这表明,回灌温度越高,越不容易造成抽水井发生热突破,因此所需的井间距也相应越小,原因在于,回灌水温度越高,同时经过热储层围岩的热补给,回灌水的温度上升就相应越快,低温回灌水对开采井的影响时间及范围也就越小。

同时,对3种不同回灌温度条件下的模拟水位数据进行了分析,如图12所示。

图12 供暖期末不同回灌温度下采灌井间剖面线水位对比Fig.12 Comparison diagram of profile line water level between production and irrigation wells at different reinjection temperatures at the end of heating period

随着回灌水温度的升高,开采井水位基本无变化,回灌井水位略有上升,但变化幅度随温度升高而减小,原因在于回灌水温度越低导致周围含水层温度降低,从而影响地下水的粘滞系数及渗透系数,回灌水温度越高,回灌井附近渗透系数越大,因此水位上升幅度也越少。

通过将渗透率及岩石热传导系数增大或减小50%进行模拟分析得出敏感度指数计算见表7。

表7 各参数不同变幅情况下敏感度指数Tab.7 Sensitivity index under different amplitude of each parameter

渗透率的敏感度指数大于岩石热传导系数,回灌井受两种参数变化的影响更大。采灌区地下水位对渗透率变化比岩石热传导系数变化更为敏感,而采灌区温度场对岩石热传导系数变化比渗透率变化更为敏感。采灌井系统对回灌水温度变化的敏感性较小,受回灌水温度变化影响仅限于回灌井附近的低温区范围及水位。

4 结语

(1)运用以TOUGH2为核心的PetraSim软件建立郑州市主城区和东部新城区地下水热耦合数值模型,取研究区水位监测点的数据以及研究区流场温度场进行模型识别验证。拟合结果较为理想,模拟结果可以反映地下水流场特征,流场及温度场与实际情况相符,可用于解释地下水流动模式。

(2)分别对研究区明化镇组及馆陶组热储层设计不同采灌方案进行模拟,结果表明在采灌方式一致时,随着井间距增大,热储层供暖期末回灌井处水位呈现上升的趋势,开采井中心处的水位降深及水位漏斗面积增大。在井间距一致时,“一抽两灌”相较于“一抽一灌”条件下,热储层供暖期末回灌井与开采井中心处的水位及漏斗面积影响呈现减小的趋势,开采井发生热突破风险也相应减小。明化镇组热储层相较于馆陶组热储层渗透率及孔隙度有所不同,在合理井间距布置方面略有差异,总之应用“一抽两灌”采灌方式比“一抽一灌”方式对地下水位及地下温度场影响小。

(3)运用局部分析法的因子变换法对模型馆陶组热储层渗透率、岩石热传导系数以及不同回灌温度进行敏感性分析,结果表明,渗透率的敏感度指数大于岩石热传导系数,其中当岩石热传导系数减小50%,开采井温度下降约0.035 ℃,岩石热传导系数增加50%时,开采井温度上升约0.02 ℃。当渗透率减小50%时,供暖期末回灌井中心水位上升1 m,开采井中心水位降低0.5 m,渗透率增加50%时,供暖期末回灌井中心水位降低0.34 m,开采井中心水位上升0.15 m。此外,采灌井系统对回灌水温度变化的敏感性较小,随着回灌温度的增加,回灌井中部低温影响范围减小,开采井水位基本无变化,回灌井水位略有上升。