地下水源热泵运行期间地下水流场和温度场的变化规律研究

2021-01-04徐连三刘红卫李永辉王富强

资源环境与工程 2020年4期

徐连三，华杉，刘红卫，柯立，李永辉，王富强

(1.湖北省地质局武汉水文地质工程地质大队,湖北武汉 430051； 2.湖北省地热能研究推广中心,湖北武汉 430051)

地下水源热泵系统将地下水作为热能载体,利用其温度常年稳定,不易受外界环境干扰的特性,夏季制冷,冬季取暖[1]。在国家宏观政策引导和地方政策支持下,地下水源热泵系统工程发展迅速,同时也暴露出一些问题,如抽—灌井布局、回灌井回灌、热贯通及系统运行对地下水环境的影响等。国内外一些学者针对如此多的问题开展了大量研究。如徐贵来等[2]通过实地对地下水源热泵系统运行期间地温变化的监测,分析研究岩土层温度变化规律,提出了提高系统能源的措施和建议；王佳乐等[3]应用有限元模拟软件模拟了岩溶含水层地下水源热泵运行对地下水温度场的影响；Park B等[4]通过试验观测发现注水井回灌导致地下水热弥散系数变化,从而影响含水层渗流场和温度场；张淑秘等[5]应用水热耦合模型模拟了不同布井模式下含水层温度变化和热贯通特征。地下水源热泵系统抽灌井布局,模拟预测地下水源热泵受水文地质条件因素影响,各区域地下水源热泵系统运行对地下水环境影响不尽相同,且运行模式不同影响亦不同。鉴于此,本研究以江汉油田供热移交改造项目——向中能源站地下水源热泵系统工程为例,利用MODFLOW软件的MT3D模块建立地下水流—热运移耦合数值模拟模型,模拟地下水源热泵运行期间地下水流场和温度场的变化规律,旨在为后期地下水源热泵系统优化设计,减轻系统运行对地下水环境的影响提供依据。

1 研究区概况

1.1 研究区水文地质条件

研究区位于湖北省潜江市广华区(见图1),距潜江市区约23 km,地貌上属典型的平原湖区地貌,地势西北部略高、东南部略低,海拔24～40 m。属亚热带湿润季风气候区,四季分明,多年日平均气温18.3 ℃。

图1 研究区地理位置图Fig.1 Geographical location map of the study area

整体属于汉江平原水文地质单元,200 m以浅主要为第四系松散岩类孔隙水,可划分为两个承压含水层组：第一层承压含水层埋深约13.0 m,厚度约为29.4 m,主要物质成分为粉砂、粉细砂,单井涌水量500～1 000 m3/d；第二层承压含水层埋深约46.0 m,厚度约为78.3 m,主要物质成分为圆砾、中砂和细砂,单井涌水量2 000～3 200 m3/d。

区内地下水位埋深随季节波动,枯水期地下水位埋深约为2.0 m,丰水期地下水位埋深0.5 m左右。地下水流向近似由北向南,根据水位标高初步测算区内天然水力坡度约为0.23‰。根据对工作区内地下水温度的监测结果显示,在150 m深度范围内,地下水取水口温度普遍在18.5～19.5 ℃。区内第四系地下水化学类型主要为HCO3-Ca·Na及HCO3-Ca·Mg,矿化度<1 g/L,属于低矿化度淡水。

1.2 地下水源热泵工程布置

区内布置有16眼抽水井和38眼回灌井,其中抽水井出水量120 m3/h,回灌井回灌量为50.53 m3/h,抽水和回灌均在第二层承压含水层中进行。抽水井和回灌井呈交叉分布,具体布局如图2所示。

图2 研究区抽灌井分布Fig.2 Distribution of pumping and filling wells in the study area

该系统工程只取暖,不制冷。运行模式分为运行期、蓄热期和恢复期,运行期时间为每年12月1日—次年2月28日,共3个月(90 d),所有抽灌井正常运行。蓄热期为每年的6月1日—9月30日,共4个月(120 d)。为了保证在取热运行期间被抽取热源后尽快得到恢复,避免次年取热受到影响,在系统停运后每年6—9月,与系统取热运行期间一致,所有抽灌井正常运行,将抽取的地下水充分利用夏季气温高的特点进行储存加热,回灌水补热温度为25 ℃。其余时间段(3、4、5、10和11月)均为自然恢复期。

2 水—热耦合数值模拟模型

2.1 水文地质概念模型

(1) 模拟平面范围。在系统分析区域地下水渗流场的基础上,参考能源站设计的抽灌井分布情况,以该区域为中心垂直和平行于水流方向各延伸1 km,以降低边界对大规模开采引起的地下水流场变化的影响,构成2 km×2 km的矩形区域作为本次模拟范围。

(2) 模拟垂向分层。根据前期勘查情况结合区域地质资料,本次模拟深度为地表至150 m深度范围内。模型垂直方向上共分为5层：第1层为人工杂填土、粉质粘土和淤泥质粉质粘土,概化为潜水含水层,厚度为13.2 m；第2层为粉砂、细砂等,概化为浅层承压含水层,层厚29.4 m；第3层为粉质粘土层,概化为第一弱透水层,厚度3.4 m；第4层为圆砾、中砂和细砂等,夹有薄层粉砂质粘土,整体概化为深层承压含水层,厚度78.3 m；第5层为粘土层,概化为第二弱透水层,厚度25.7 m。根据地表高程数据及各层厚度,对模型各层的顶底板高程进行赋值,构建研究区三维结构模型。

(3) 边界条件界定。区内地下水整体上由北向南径流,由此平行于地下水流向的东西两侧边界条件取为隔水边界,垂直于水流方向的北侧边界为侧向补给边界,南部边界为侧向排泄边界,两者均取为流量边界。顶部边界为开放边界,接受降雨入渗补给、农田灌溉水入渗补给、人工开采、蒸发排泄等,概化为流量边界。根据江汉平原地下水均衡计算结果可知,地下水补给量与排泄量相对均衡差为2%左右；江汉平原区地表水丰富,地表水与地下水之间交互频繁,降水量及其入渗补给量均比较大,区内地下水系统整体上保持均衡状态。因此将研究区内第四系含水系统视为动态平衡系统,在模型中将入渗补给量等总补给量与地下水开采量和蒸发量等总排泄量视为均衡,而且模型的范围比抽水回灌区大得多,还存在回灌补给,开采对边界流量的影响可忽略,仍然保持原来的进出平衡,因此在模型中取为零流量边界。模型底边界为粉质粘土层,概化为隔水边界。

研究区范围较小,大地热流的空间变异小,不同位置的地温梯度变化小,地下水以水平运动为主,水流对地温的影响较小；回灌情况下地下水开采对地下水流场的影响范围较小,对边界的影响可以忽略,因此可将侧边界取为定温度边界。上边界为该区的恒温带,设置为定温度边界,近似取值为当地多年平均气温。同样的方式,下边界亦取为定温度热流边界。

2.2 数学模型

2.2.1地下水流动数学模型

地下水—热耦合模型涉及地下水流动和热运移过程。根据渗流连续性方程和达西定律,地源热泵系统地下水流动可用偏微分方程及其定解条件表示：

(1)

式中：Kx、Ky和Kz分别为x、y、z方向的渗透系数,m/d；H为含水层的水头值,m；M为含水层厚度,m；ε为源汇项,m/d；S为给水度或比弹性释水系数,潜水含水层取重力给水度,承压含水层取弹性释水系数；Ω为模拟范围；n为边界外法向方向单位向量；Γ为侧边界；B为底边界；H0为初始水头,m。

2.2.2热量运移模型

地下热水主要以液态的形式存在,在多孔介质中地下热水热量运移模式主要包括热跟随水流传输(对流)及热量在水和多孔介质中的热传导,其方程可表示为[6]：

(2)

式中：n表示孔隙度；ρm表示含水介质的密度,kg/m3；ρw表示液体的密度,kg/m3；α表示纵向弥散度,m；cw表示液体的比热容,J/(m3·K)；cm表示含水介质的比热容,J/(m3·K)；λm表示含水介质热传导系数,J/(m·s·K)；qh表示热量的输出或输入,J/s。

由于地下热水系统中的热量传输机制以对流、传导为主,其数学表达式与地下水系统中溶质运移的数学模型一致。而且地下水溶质运移的数值模拟程序更为成熟,有学者使用地下水溶质运移模拟程序开展地下热水的运移和热传输模拟,模拟结果显示MT3DMS的模拟结果与FEFLOW模拟结果相近。通用地下水数值模拟软件(GMS v10.1)教程中也包含使用MT3DMS进行地下热水流和热传输模拟。

溶质运移的数学方程为：

(3)

式中：ρb表示含水介质的密度，kg/m3;Kd表示有效分布系数,m3/kg；Dm表示有效分子扩散系数,m2/d；va表示地下水流速,m/d；qCk表示源汇项,kg/(m3·s)。

从方程式(2)和(3)可以看出,MT3D模块求解溶质运移与热量运移方程的形式是一致的,MT3D模块中溶质运移与热量运移参数之间的关系如下式所示：

(4)

(5)

(6)

2.3 数值模型及参数设置

2.3.1离散化处理

在水文地质概念模型基础上,利用MODFLOW数值模拟软件的MT3D模块建立数值模拟模型。本次模拟范围为2 000 m×2 000 m×150 m。为了提高模型的计算速度,同时保证模型的计算精度,平面上采用有限差分法进行等间距剖分,网格间距为50 m。水平方向上,在抽水井和回灌井群区范围内,对剖分网格进行加密处理,网格间距为10 m。模型垂直方向上共分为潜水含水层、浅层承压含水层、第一弱透水层、深层承压含水层和第二弱透水层等5层,共56 160个单元。

本次模拟以能源站热源井刚开始运行的时间(2019年12月1日)作为初始时刻,以天为单位,对能源站地源热泵系统运行20年后的地下水流场和温度场特征进行预测。地下水三维水—热模型计算过程中均采用非稳定流模拟,模型计算初始步长设置为1 d。

2.3.2边界条件

(1) 水流边界。根据研究区能源站钻孔抽水试验计算结果可知,抽水量为50 m3/h、80 m3/h和130 m3/h时,计算的影响半径R分别为56、100和300 m。研究区侧向边界范围远远大于抽水井影响半径,因此将侧向边界作为隔水边界进行处理,取零流量。模型底边界为隔水边界,取零流量。

(2) 温度边界。将模型第一层设置为恒温带,温度取潜江市的多年日平均气温,取值18.3 ℃。

侧边界取为定温度边界,边界温度根据地温梯度进行计算,赋入模型中。参考《中国大陆地区大地热流数据汇编(第四版)》中江汉盆地的地温梯度,取值范围19.7～27.9 ℃/km,平均值为22.9 ℃/km。底边界取定温度边界,温度按照地温梯度计算而来。

2.3.3初始条件

本次模拟主要根据研究区的地表高程数据,结合向中能源站钻孔时监测的水位埋深,对模型中各剖分单元的初始水位进行赋值。研究区热储层(深层承压含水层)地下水初始水位高程取值30 m。

依据地表多年日平均气温数据,结合研究区地温梯度,按照各层的厚度计算其温度,然后作为初始温度分层赋入模型。研究区热储层(深层承压含水层)初始地下水温度取值19.5 ℃。

2.3.4模型参数的选用

构建三维水—热耦合运移模型时,需要对水文地质参数和含水层骨架的热物理性质参数进行赋值。根据前述水文地质概念模型设计,本次模拟需要对渗透系数(Kxx=Kyy、Kzz)、给水度、弹性释水系数、孔隙度、流体热传导系数、固体热传导系数、流体比热容和固体比热容等参数进行赋值。根据野外试验和室内热物理参数试验,结合相关文献资料[7],初步确定了模型所需参数,不同层位(潜水含水层、浅层承压含水层、第一弱透水层、深层承压含水层和第二弱透水层)参数数值如表1所示。

表1 模拟区各层参数统计表Table 1 Statistical table of parameters of each layer in simulation area

利用MT3D模块进行溶质运移和热量运移模拟时,主要利用Advection、Dispersion和Chemical reaction程序包进行方程求解。其中Advection程序包采用软件默认的参数进行,Dispersion和Chemical reaction程序包中各层参数赋值如表2所示。

表2 模拟区Dispersion和Chemical reaction程序包中各层参数表Table 2 Parameter table of each layer in the simulation area Dispersion and Chemical reaction package

2.4 抽灌井运行模拟结果分析

在验证模型边界条件、初始条件、参数分区和源汇项的基础上,按区内地下水源热泵工程实际运行模式,建立预测模型,模拟热源井运行20年的地下水流场和温度场的空间分布。

(1) 地下水流场特征。地源热泵系统供暖期地下水流场分布如图3所示,热源井进行供暖期间在抽水井分布比较集中的区域形成小型水位降落漏斗,与初始地下水位(30 m)相比,此时水位落差接近4.5 m。实际抽水和回灌模式下热源井进行完全回灌,水源热泵系统运行1年和20年后的流场空间分布趋势相同,并均与研究区初始地下水流场空间分布保持一致。

图3 热源井运行90天的地下水位等值线图(图中红色圆点代表抽水井,黑色圆点代表回灌井)Fig.3 Contour map of groundwater level in 90 days operation of heat source well

从图4可以看出,地源热泵系统刚开始运行时,抽水井处水位迅速下降,此后水位保持较稳定状态,停止运行后,水位迅速恢复至初始状态,随着夏季蓄热期地热井的运行,抽水井处水位同样迅速下降后保持稳定,恢复期水位迅速恢复。回灌井处的水位随着各个时期的变化呈现先升高、在保持稳定后下降的趋势,如此循环反复(图5)。

图4 不同位置水位动态随时间变化曲线(运行1年)Fig.4 Water level dynamic curve with time at different positions

图5 不同位置水位动态随时间变化曲线(运行20年)Fig.5 Water level dynamic curve with time at different positions

(2) 温度场特征。地源热泵系统供暖期的地下水温度场分布如图6所示。地源热泵系统在运行期间,通过井群回灌向地下含水层注入冷水,此时注入的冷水堆积形成以回灌井为中心的冷水体,然后向回灌井周边延伸,导致回灌井周边温度降低。系统运行期间距离回灌井10 m和20 m处地下水温度受回灌冷水体的影响呈现迅速下降的趋势,温度分别降低至9 ℃和12 ℃；夏季蓄热期(6—9月)回灌水体温度较高(25 ℃),此时温度呈现迅速上升的趋势,温度升高至24 ℃和22 ℃,恢复期温度维持在较稳定状态。距离回灌井40 m处和抽水井附近受回灌冷水体影响相对较小,温度呈现缓慢下降趋势,与初始温度(19.5 ℃)相比,温度下降2 ℃左右；夏季蓄热期，地下水继续受回灌冷水体影响呈现下降趋势,之后恢复期受蓄热期回灌水体的影响,温度呈上升趋势,并升高至18 ℃左右。

图6 热源井运行90天的温度等值线图Fig.6 Temperature contour map of heat source well in 90 days operation

抽水井和回灌井呈交叉分布,而且存在夏季蓄热期(6—9月),此时地源热泵系统运行1年和20年的温度场如图7和图8所示。结合不同位置温度随时间变化动态曲线(图9和图10)可以看出,该模式下地源热泵系统在回补(地下热蓄)期末抽水井处温度回升至接近20 ℃,回灌井处温度稳定在24 ℃左右。由此可以看出,当抽水井和回灌井交叉分布,且存在蓄热期时,抽水井和回灌井处温度可达20～24 ℃,这种模式对提高系统能效和保证冬季供暖效果更为有利。