潮白河春季生态补水及地下水响应

2022-04-25李志强王亚娟温子希张秋兰崔亚莉张琳琳杨国军张莉娜刘翠珠董鹤

南水北调与水利科技 2022年2期

李志强,王亚娟,温子希,张秋兰,崔亚莉,张琳琳,杨国军,张莉娜,刘翠珠,董鹤

(1.中国地质大学(北京)水资源与环境学院，北京 100083;2.北京市水文总站,北京 100089)

地下水资源是水资源的重要组成部分。对于北京而言，地下水资源既是首都发展的基础性资源，更是北京市发展的战略资源和重要保障[1-4]。潮白河冲洪积扇中上游地区是北京市地下水资源最丰富的地区，在城市供水中占据了重要地位。由于多年连续超采,地下水位持续下降[5-8]。2014年12月南水北调中线工程正式通水，截至2021年3月已累计向北京供水超过62亿m3，在缓解北京市供水紧张局面的同时，也为地下水的涵养提供了契机[9-11]，2015年通过小中河向怀河及潮白河牛栏山地区补水，通过北上台及大水峪水库对潮白河地区进行试验性补水，共引水3 700多万m3；2016—2019年均通过南水北调中线工程引水及水库放水对潮白河流域进行补水。为贯彻落实京津冀协同发展战略，深入落实北京城市总体规划，大力推进北京市委市政府“用生态的方法解决生态的问题”、积极推进“藏水于地”的部署要求，2021年春季对潮白河实施了1个月试验性生态补水，本文即基于此次生态补水展开相关研究。

南水北调水进京后，在满足居民生活基本用水的情况下，按照“喝、存、补”的原则，部分来水通过水库放水、自然河道补水的方式对水源地进行回补。补水后地下水水位、水量、水质各方面对生态补水动态响应的研究成为关注重点。比如：永定河生态补水工程中，李海军等[12]和汤万龙等[13]利用监测数据分析补水期间水位、水质的变化，运用水均衡法分析生态补水对地下水储存量的影响；胡立堂等[14]对永定河补水后地下水动态变化和补水滞后性进行分析，划分了观测孔地下水动态类型，为补水方案提供参考；霍丽涛等[15]则以2015—2016年潮白河调水回补地下水为例，分析了观测孔水位变化和水质变化。除了利用监测数据分析补水对地下水的影响外，学者们还通过建立模型的手段对补水动态响应进行研究，其中：朱冰等[16]、彭勃等[17]和刘波等[18]采用Visual Modflow软件构建黄河三角洲地下水流数值模型，模拟分析湿地生态恢复区生态补水对地下水的影响；Sun等[19]基于人工神经网络等机器学习建立地下水模型预测了永定河补水入渗条件地下水位变化。可见构建地下水流模型能较好地实现对地下水动态的模拟，有利于分析生态补水效果。GMS是由多种模块构成的一款可视化地下水数值模拟软件，处理功能强大，适用范围广泛，是使用最多的地下水模拟平台之一[20-21]。本文采用GMS软件构建潮白河流域地下水流数值模型，模拟分析地下水对2021年春季潮白河生态补水的动态响应，以期为科学实施生态补水工程提供参考。

1 研究区概况

研究区位于北京市东北部，包括密云、怀柔、顺义3个区县的部分平原区，研究范围及边界条件见图1，总体地形西北高、东南低，地面由北向南倾斜。研究区属温带大陆性半湿润季风气候，多年平均气温为11.7 ℃，日极端最高气温可达41.3 ℃，日极端最低气温可低至-23.3 ℃。年均(2000—2020年)降水量为542.58 mm，见图2，全年降水最多的月份为6—9月，降水量占年降水总量的82%以上。

图1 研究区范围及边界条件

图2 2000年至2020年苏庄站降雨量

研究区内主要水系为潮白河水系。上游为潮河和白河。白河发源于河北省沽源县，流经赤城县，进入北京境内，由延庆、怀柔汇入密云水库；潮河发源于河北省丰宁县，经栾平、密云注入密云水库。潮河、白河出库后在密云县河槽村汇合为潮白河[22]。根据前人研究成果[23]，依据地貌单元、地层沉积环境、地下水的赋存条件以及水力联系等特征，研究区属于潮白河地下水子系统，在平原区以及山间沟谷地带普遍分布有第四系松散沉积物，主要特征为北薄南厚，东薄西厚；颗粒从北向南由粗变细；层次也由单一到多层。根据地下水埋藏条件、水动力特征和水质特征，以及该区水文地质和开采利用条件，把第四系松散沉积岩类细分为单一结构的砂卵砾石层、2～3层结构的砂卵砾石层、多层结构的砂砾石层夹少量砂层。

2 数据及研究方法

2.1 数据基础

地下水模型构建主要需要钻孔资料、补给项、排泄项和模型校准数据等。为构建长时间序列地下水模型，主要数据如下。

收集研究区DEM高程数据、100个钻孔资料和5条剖面，结合开采层位选取相对稳定厚度较大的隔水层顶板作为分层界线来确定主要含水层的顶底板标高。

含水介质参数和水文地质参数初值根据抽水试验、相关水文地质报告及前人研究成果[23]确定。

密云、怀柔、顺义3区县的气象站2000—2020年日尺度降雨量数据。

密云、怀柔、顺义3区县2000—2020年工业、生活、农业年开采量；怀柔应急水源地、水源八厂、潮-怀应急水源地等集中水源地开采量。

密云、怀柔、顺义3区县内长观孔位置及2000—2020年逐月水位埋深。补水期间(2021年4月30日—5月27日)监测井日尺度地下水位监测数据。

南水北调2015—2021年补水路径及补水具体情况数据。

数据主要来源为北京市水文总站，其他未说明数据均采用前人研究数据[23]。

2.2 研究方法

2.2.1模型原理

GMS(groundwater modeling simulation)软件是国内外常用的地下水流模拟软件，其包含的MODFLOW程序，是基于有限差分法刻画地下水流运动规律的计算机程序。此程序具体求解方法是通过把研究区域进行空间和时间上的离散，建立每个网格单元的水均衡方程式，所有网格单元的方程联立成大型的线性方程组，耦合初始条件及边界条件后，迭代求解方程组即可求得每个网格单元的水头值。GMS 10.0版本以后，嵌入了MODFLOW-USG(MODFLOW-Unstructured Grid)程序，该程序基于体积有限差分法对网格进行剖分[24-25]。区别于以往剖分方式，MODFLOW-USG程序对单元格的定义不再局限于行、列、层，剖分方式更加灵活多变，可以使用任何类型的网格，即剖分时具有任意的面和节点的单元格，这种特性可以使得尖灭等地质特征进行更真实建模，这也促使在进行网格剖分时可以对重点关注区域(如河道、集中水源地)进行局部加密处理。

2.2.2模型建立

根据对研究区水文地质等条件的分析，确定模拟范围为北起密云水库山前一线，南至顺义区苏庄；东西界限为顺义行政区划界限各外推约1 km2。模型总面积约1 440 km2。

研究区北、西北、东北3个方向接受山区侧向补给，为自然边界，边界流入量根据资料确定；南部与平原区接壤，与平原区的地下水发生水量交换，为通用水头边界，边界流入流出量根据模型通用水头包计算得到。潜水面为该模型的上边界，地下水通过此边界得到降水入渗补给、山前侧向补给、河水渗漏补给、灌溉回归入渗补给、渠系管网等入渗补给、生态补水等。模型的底边界为基岩，故将其看作不透水边界。

此次模型网格加密方法选择的四叉树网格剖分，对补水河道、重点关注的水源八厂、怀柔应急水源地等地进行了加密。其中,基础网格为200 m×200 m，细化网格100 m×100 m。空间上将模型分为四层:第一层为潜水含水层;其他层为承压水含水层。模型期为2000年1月至2020年12月，以1个月为1个应力期，共252个应力期；选择2000年1月至2011年12月为模型识别期，2012年1月至2020年12月作为模型验证期，以2000年1月流场作为初始流场，建立长时间序列地下水流模型。

模型中面状补给项包括降雨入渗补给和灌溉回归渗漏补给，面状排泄项包括蒸发和部分人工开采(如工业、生活、农业开采)。降雨入渗补给根据每个降雨入渗分区把降雨量处理成为降雨强度，降雨入渗分区依据降雨入渗系数不同划分，利用RCH(recharge)程序包处理输入模型。蒸发项则将蒸发量处理为蒸发强度，用EVT(evapotranspiration)程序包处理输入模型。渗漏补给和部分人工开采则以注水井或抽水井的形式将数据处理成相应格式用井(Well)程序包输入模型。线状的补排项主要有山前侧向补给量、南水北调补给水量和部分通用水头边界流入流出量，如山前侧向补给量主要应用达西(Darcy)定律计算得到每个网格的补给量，应用ArcGIS将数据拓扑到相应位置网格中心点上，并计算转换成单个网格上的井流量，利用Well模块处理。通用水头边界的流入流出量根据Darcy定律和流场形态，通过GHB(general head boundary)程序包计算得出后代入模型。

2.2.3识别及验证

利用同时期的地下水流场或者长观孔历时曲线来拟合模拟水位，识别地下水流模型中水文地质参数、边界值和其他源汇项，使模型愈加符合研究区的水文地质条件，以便更精确地定量研究模拟区的补给与排泄，预报给定水资源开发利用方案下的地下水位变化。利用区内4个长时间序列观测孔对地下水流模型进行参数率定，模型主要参数率定结果见表1。选取国内外研究常用的评价因子决定系数(R2)和均方根误差(RMSE)作为评价模型参数优化的评价标准。R2越接近1，表明模拟值越接近实测值，模拟效果越好;RMSE则反映了模型的模拟值与真实观测值的接近程度[26]。模型验证结果(图3)显示：4个观测孔决定系数范围为0.93～0.97；NH1、NH2、NH3三孔的RMSE范围为3.19～3.78，另一个孔的RMSE为6.97，受开采和边界流入量影响，数值偏大。评价因子R2和RMSE值表明模型源汇项、参数基本合理，验证了模型的可靠性，模型可以用于生态补水过程分析。

图3 长时间序列模型水位拟合结果

表1 模型主要参数率定结果

3 潮白河流域生态补水的模型应用

3.1 潮白河春季生态补水概况

潮白河流域试验性补水在2014年底开始实施，利用“密云水库调蓄工程”，在汛期通过小中河向怀河及潮白河牛栏山地区补水，通过北上台水库和大水峪水库向潮白河地区补水。2018年9月，南水引水及水库放水累计引水量达到4.56亿m3；截至2019年9月，南水引水及水库放水累计引水量超过5.7亿m3。2021年4月30日起，潮白河流域实施22年来的首次较大规模试验性补水，主要通过白河、潮河、雁栖河、怀河、小中河向白河干流集中补水，计划补水总量达到2.4亿m3。截至2021年5月27日，累计出库及南水北调放水总量达2.06亿m3。

3.2 潮白河春季生态补水模拟

在长时间序列地下水流模型的基础上，将模型应用于潮白河流域2021年春季试验性生态补水中，模型边界条件及参数不变，预测期为2021年4月30日—2021年7月31日，1 d为1个应力期，共93个应力期；初始流场为2021年4月末流场。补水期内(2021年4月30日—5月27日)无降水，为了充分体现地下水对生态补水的响应过程，补水期结束后的模型应力期内也没有考虑降水的影响，故模型补给项主要包括灌溉入渗、侧向流入、水库渗漏及潮白河对河道补水；研究区内长期超采至今，导致地下水埋深大于蒸发极限埋深(4 m)，故排泄项未考虑蒸发，主要为人工开采。将数据处理成相应的格式输入到GMS中，利用模型分析评价潮白河2021年春季生态补水效果和地下水响应。

4 结果与讨论

4.1 地下水水位响应

本次补水沿潮白河流域10 km范围内共布置了65眼观测井作为补水监测井，对所搜集的监测数据进行水位拟合，误差基本在0～1 m的可接受范围内，其中4眼观测井拟合情况见图4。模型模拟水位和实际水位数据较为接近，模型参数和源汇项基本符合实际情况。分析拟合孔水位变化可知：7号孔反应迅速，补水后不久水位即开始上升，这与该孔位于上游且离河道较近有关；其他拟合孔水位在补水一段时间后，水位上升明显。截取剖面A-A′，见图5，随着时间的推移，到7月31日距离主河道两岸的水位开始消散，而距离河道较远的地方水位逐渐回升。24号和27号拟合水位出现0～1 m波动，拟合效果一般，可能是在补水初期受到怀柔应急水源地和水源八厂开采的影响。

图4 观测孔水位响应及拟合结果

图5 A-A′剖面水位变化

4.2 补水影响范围

根据模型运行结果，考虑生态补水作为单一变量影响下地下水水位动态的变化。在补水预测模型基础上，考虑研究区无补水情形下的地下水水位。将两者补水期末的水位相减得到补水期间动态水位抬升图，见图6。分析可知，补水期末，沿着补水路径周围水位均有所上升，最大上升幅度接近18 m。以水位抬升值大于0划定此次补水的影响范围，计算得到研究区内的影响面积约为842.8 km2，占研究区总面积的58.6%。靠近潮白河干流沿岸地下水位上升幅度大于10 m，影响面积约为44.2 km2，占比约为3%。

图6 补水期末水位变幅

4.3 补水对储存量的影响

补水期间，补水对地下水的入渗是研究区的主要补给来源。根据补水期间地下水观测资料，确定观测孔初始水位值和补水期末水位差值，圈定各观测孔所控制的面积，取用相应的给水度，算出整个研究区的储存量变化值。

(1)

式中：ΔQ储为储存量的变化量，m3；μi为第i个观测孔所对应的给水度；Δhi储为第i个观测孔的水位差值，m；Fi为第i个观测孔所圈定的面积大小，km2。

根据式1计算可知,在总水量为2.06亿m3的前提下，补水渗漏补给地下水量为1.57亿m3(约占总水量的76%)，研究区第一含水层系统储存量增加了1.49亿m3。由此可见，生态补水对于潮白河冲洪积扇地区地下水的恢复具有重要的意义。

4.4 地下水位的滞后效应

地下水的抬升和生态补水有着密切的联系，但是地下水水位的回升需要一定的时间，即存在一定的滞后性[26]。根据现有监测井分布特点，选取距离潮白河河道两岸的4个观测井分析地下水水位对补水的滞后响应，见图7,南岸为13号(3 km)、17号(6 km)，北岸为11号(3 km)、19号(6 km)。对比同侧观测井可以发现：13号、11号观测井离河道较近，对补水的响应及时，水位抬升迅速，响应滞后时间较快(t0=3 d,t2=6 d);17号、9号距离河道较远，地下水响应时间较慢(t1=30 d,t3=18 d)，水位上升趋势较平缓。13号观测井响应滞后时间最短，水位上升最多，分析其原因是从补水初期开始最先受到补水入渗持续的影响，17号观测井响应滞后时间最长，水位上升缓慢，其原因是距河道较远且仅受到潮白河干流补水的影响。