鄱阳湖典型洪泛区地下水数值模拟研究*
2023-01-13曹思佳李云良李宁宁宋炎炎赵贵章李志萍
曹思佳,李云良,李宁宁,陈 静,宋炎炎, 赵贵章,李志萍
(1:华北水利水电大学地球科学与工程学院,郑州 450045) (2:中国科学院南京地理与湖泊研究所,南京 210008) (3:河海大学水文水资源与水利工程科学国家重点实验室,南京 210009) (4:江西省地质局工程地质大队,南昌 330000) (5:江西省水文监测中心,南昌 330002)
洪泛湿地是河流与陆地之间典型的水文过渡带,也是一种受洪水干扰频繁的特殊下垫面类型,约占全球湿地面积15%[1]。从水循环角度而言,洪泛湿地是大气降水、地表径流以及地下水相互转化的重要界面,因洪水过程的干扰,地表-地下水之间相互转化通常具有高度的动态性和敏感性[2]。因此,洪泛湿地系统具有水文和生态环境研究的多重特色和意义,地表和地下水文过程的季节性变化对水流界面交互、营养物质交换以及沉积作用有重要影响作用,很大程度上促进了洪泛系统物质流、能量流和信息流的转化和传递[3]。诸多研究表明,河流或者湖泊的水文情势变化以及地表土地利用类型等对洪泛湿地地下水具有重要影响和贡献,主要体现在地下水位的动态性和水动力场的变异性以及地下水对湿地的动态补给模式等[4]。地表和地下水文研究已成为全球洪泛湿地研究的重要问题,地表-地下水交互转化是生态系统研究的基础,是研究水环境防治的重要组成部分,同时也对洪泛湿地的管理和保护起至关重要的作用[5]。
在我国的长江中下游地区,以湖泊为核心的多类型复合湿地生态系统分布颇为广泛。其中,鄱阳湖及其洪泛湿地系统具有水文和生态环境等诸多领域研究的区位优势,也是国际上公认的重要湿地。鄱阳湖洪泛湿地在维护系统平衡、植被群落分布等方面发挥着不可替代的作用和意义,也是我国重要的候鸟越冬栖息地[6]。然而,由于近年来人类活动的干预以及气候条件变化的叠加影响,鄱阳湖水位呈现出洪季偏洪、枯季偏枯的总体态势,尤其是鄱阳湖低枯水位的屡次出现,无疑干扰了湖泊洪泛的自然水文条件和生态状况,导致洪泛湿地的淹没动态、出露面积和地表形态等诸多因子发生显著变化,将会使鄱阳湖洪泛湿地面临外部复杂环境条件及水资源转化的不确定性[7]。因而,以湖泊变化环境条件为背景,聚焦湖泊洪泛湿地地表-地下水转化及其水量动态研究,对区域水资源优化配置、湿地系统演变机制以及湿地生态系统健康等具有重要理论意义和实际应用价值。
近年来已有学者针对鄱阳湖洪泛湿地系统的地表和地下水方面开展了一些探索性研究工作,取得了一些重要进展和认识。先前研究发现,鄱阳湖水位变化对洪泛湿地系统的地下水动力场和湖水-地下水转化关系具有不可忽视的影响作用,同时也对洪泛地表生态植被等方面具有重要影响[8]。Li等[9]基于野外观测数据和水化学分析,发现了鄱阳湖与洪泛湿地的地下水之间存在密切转化关系,但受资料和方法所限,还无法深入认识地下水对湖水变化的动力响应过程以及两者之间的转化通量。陈静等[10]联合应用野外观测资料和能量守恒等计算方法,估算了鄱阳湖洪泛区典型碟形湖与地下水的交换通量,但由于原位监测无法精确刻画研究区的空间变异情况以及资料连续性相对较差,通量结果的估算误差还无法确定。Zhang等[11]利用氢氧稳定同位素示踪技术分析了鄱阳湖洪泛湿地的土壤-植物-大气连续体(SPAC)的水分传输,但研究结论仅是定性得出地下水对湿地水流迁移具有一定的促进作用,尚无法定量识别地下水与地表水之间的转化与动态交互。Wang等[12]应用地下水数值模型开展了鄱阳湖赣抚平原区的地表-地下水相互作用研究,分析了地下水对河湖水体的贡献作用及对气候变化和人类活动的响应,但相关工作主要侧重于流域下游的赣抚平原地区,对洪泛区地下水动力学过程和转化机制的认识相对不足。上述分析可知,目前围绕鄱阳湖洪泛系统地表-地下水转化等方面已取得了一些共性认识,强调了地下水文的动态性及对湿地生态的贡献,但研究方法以原位观测和示踪技术为主,开展洪泛区地下水动力场及其响应过程的数值模拟研究对深入理解和完整诠释湿地生态环境效应内涵具有重要意义。
在当前变化环境影响下,鄱阳湖洪泛湿地的水文和生态环境问题愈发突出。本文以鄱阳湖典型洪泛湿地为研究区,基于野外观测数据和资料获取,应用地下水数值模型,着重从地下水对湖水的贡献作用的视角来分析鄱阳湖洪泛区的水资源转化问题。本文主要研究目标为:(1)构建基于FEFLOW的洪泛湿地二维地下水动力学模型,分析地下水位的时空动态变化与特征;(2)定量识别湖泊不同水文时期地下水动力场的演变规律及对湖水位变化的响应过程;(3)基于地下水均衡分析,揭示不同水源组分对洪泛区地下水系统的相对贡献。本研究可为后续水文地球化学循环、生态环境效应等相关研究提供重要数据资料和科学支撑。
1 材料与方法
1.1 研究区
鄱阳湖位于江西省北部及长江中下游南岸,是一个过水性、吞吐型湖泊,主要承接赣江、抚河、信江、饶河和修水五河来水,经调蓄后由北部湖口注入长江干流(图1a)。鄱阳湖流域南北长约173 km,东西最大宽度约74 km,北部入江水道最窄处宽度约2.8 km,湖岸线总长达1200 km(图1c)。鄱阳湖流域面积占长江流域面积的9%,平均年径流量为1450亿m2,水量约占长江流域水量的5%[13]。鄱阳湖流域三面环山,周围高中间低,地形高程变化为30~2200 m。鄱阳湖地处亚热带暖湿季风气候区,年内降水量分布不均,加上流域五河来水变化,导致其具有高度变异的湖泊水位变化特征[14]。湖泊洪、枯季节的水位差可达15 m,形成“洪水一片,枯水一线”的独特洪泛特征。洪水期间湖区整体被淹没,平均水深不足6 m,最大水深可达30 m,湖水面积可达3000 km2,而枯水期水域面积则萎缩至不足1000 km2[15],大部分水体基本萎缩在主河道附近。在鄱阳湖洪泛湿地内部,由于水文水动力、泥沙冲淤等原因形成众多大大小小、形状各异的洼地,共计约102个,这些洼地被称为碟形湖[16]。在洪水期,这些碟形湖被洪水淹没,与主湖融为一体,而在枯水期,湖水位快速下降,这些碟形湖逐渐形成相互独立的子湖泊水体。碟形湖为洪泛湿地生态环境提供了重要水和物质保障。碟形湖附近的植被长势良好,碟形湖草洲面积占整个湿地草洲面积的23.09%[17]。同时碟形湖中生存着大量的底栖动物,为候鸟越冬提供了物质基础和适宜的生存条件。鄱阳湖受动态水情波动的影响,主湖区和洪泛湿地并没有明确的界限,但普遍认为常年有水的地方为主湖区,受主湖区水位变化影响的周边广泛滩地则视为洪泛区。根据鄱阳湖洪枯季节的水面积变化范围,本文典型洪泛区主要位于鄱阳湖主河道(星子-都昌-康山)和西侧湖区边界之间的洲滩湿地,面积约为1646 km2(图1c)。
图1 鄱阳湖流域主要水系分布(a)、湖区周边流域地层岩性(b)和典型洪泛区观测站点分布(c)Fig.1 Major rivers in the Lake Poyang catchment (a), formation lithology around the lake (b) and hydrological gauging stations distributed in the typical floodplain area of the lake (c)
1.2 基础数据与获取
鄱阳湖湖盆地形高程数据原始分辨率为5 m×5 m,主要用来刻画洪泛区地形特点及构建地表-地下水数值模型。星子、都昌和康山水文站的日水位观测数据来源于江西省水文局,主要用来描述湖泊水位的季节性变化及作为数值模型的边界输入条件。地下水观测井(1#~5#)的日水位数据来源于加拿大Solinst传感器的野外自动观测,观测点主要分布在吴城国家级自然保护区湿地以及南矶保护区湿地,该数据主要用来反映地下水位的波动状况及验证地下水模型。在修水和赣江的下游平原地区(未控区间),鄱阳湖湖区主要岩性为第四纪松散岩类沉积物,周围以变质岩为主,零星分布岩浆岩、碳酸岩及少量碎屑岩(图1b)。未控区与湖区西侧的交换通量根据未控区的地下水位资料和达西定律计算获取,未控区地下水位数据来源于江西省水文局,通量数据用于作为模型的边界输入条件。降雨和蒸发的日观测数据来源于中国科学院南京地理与湖泊研究所鄱阳湖湖泊湿地综合研究站(星子站附近),用于作为数值模型的大气输入条件。典型洪泛区碟形湖与地下水之间的日通量变化数据采用先前的估算结果[10],作为数值模型的地表源汇条件。文中涉及的一些水文资料的具体观测点位及分布可参考图1c。考虑到鄱阳湖水情变化特征的普适性以及所需不同数据资料的完整性,故选取2018年作为研究期来开展正常年份条件下的地表-地下水文过程分析。
1.3 数值模型介绍
FEFLOW(Finite Element subsurface FLOW system)最初是由德国水资源规划与系统研究所WASY公司开发的地下水有限元数值模型,目前已成为丹麦水利研究所DHI模型系统的重要组成部分。模型采用有限单元离散技术,应用适应性更好的非结构性网络,以此实现局部网格加密及其重点区域的精细化刻画,进而来灵活应对各种地下水流场模拟与溶质运移模拟中可能出现的复杂的物理过程[18]。FEFLOW可以实现多孔介质达西流、非饱和流、潜水水流模拟和迁移、变密度流和裂隙等诸多实际问题。因此FEFLOW模型已被广泛应用于地下水动态预测、地下水资源利用分配、水热耦合运移、地下水污染运移、由于抽水引发的地面沉降以及海水入侵方面的研究,并取得了许多重要成果[19-22]。其地下水流运动数学模型为:
(1)
1.4 地下水模型构建
根据鄱阳湖洪泛区DEM地形高程特征,本文共计提取出41124个高程点,采用克里金插值方法,将插值结果作为数值模型的实际地表高程。结合湖区地形地貌分布格局,研究区共计剖分6329个三角形有限单元网格和8680个节点,三角形边长变化范围介于20~2000 m之间,很好地刻画了洪泛区的复杂地形特点(图2a)。考虑到本文典型洪泛区实属鄱阳湖区的一部分,相对大尺度流域而言,其地质类型和成因相对一致,加上洪泛区野外资料极其有限,其中一些水文地质参数未考虑分区特征(比如降雨入渗系数、渗透系数和给水度等)。根据研究区周边的野外钻孔资料(图2b),湖岸带岩性主要为细砂和粉质黏土,粉质黏土中的粘粒含量较少,与细砂的岩性特征相差不大,且根据现场抽水试验结果,细砂和粉质黏土层的渗透系数较为接近,在地表以下20 m处有一黏土层,透水性较差,其埋深情况相对均匀,可作为含水层底板。因此本文将地下含水层概化成一层,将潜水含水层底板埋深设为地表以下20 m处。研究区的东侧边界是鄱阳湖主河道,从北部星子延伸至南部康山,主要受湖泊水位变化控制,因此模型东侧设置为给定水头边界条件,根据星子、都昌和康山水位观测资料进行插值并分为3段来分别给定,以体现东部边界的湖水位空间差异(图2a)。洪泛区湿地的西侧主要接受修水和赣江两大区间的地下水补给作用,因此西侧边界根据修水和赣江的影响范围,划分为两段并分别设置为给定流量边界(图2a)。
图2 研究区模型概化示意图(a)和都昌湖岸钻孔柱状图(b)(钻孔数据来源于江西省地质局)Fig.2 Conceptual map of groundwater flow model (a) and vertical borehole profile near the lakeshore of the Duchang gauging station (b) (data obtained from Jiangxi Geological Bureau)
研究区主要的源汇项为大气降水、地表蒸发及碟形湖入渗补给或排泄。如前所述,本文选取2018平水年作为模拟期,降水蒸发等气象条件与多年平均情况相当,具有很好的代表性。根据先前文献资料,研究区的入渗系数取值范围为0.06~0.1[23],空间变化总体上较小,故本文将研究区的降雨入渗补给系数给定为0.1,即认为10%的有效降雨入渗补给潜水含水层。参照许秀丽等[24]的已有研究,鄱阳湖洪泛区湿地的土壤年蒸发量大概为100~200 mm,约为E601B型蒸发皿年蒸发总量的10%,考虑到植物蒸腾作用主要发生在每年3-5月且主要吸收土壤水,地下水的贡献比重相对较小。因此,本文以E601B型蒸发皿数据为基础,乘以折算系数0.1作为地下水的蒸发量。相对于鄱阳湖周边的地层岩性分布而言,本研究所覆盖的洪泛区岩性分布整体上较为均一(图1b),加上研究尺度相对较小,因此洪泛区的地下水数值模拟参数设置为均质各向同性。根据2018年鄱阳湖平水年遥感影像图(Landsat),对数值模型计算域内的52个碟形湖进行定义和设置(图2a)。在野外现场调研中,发现碟形湖底部存在明显的淤泥质弱透水层,其厚度约为0.3~0.6 m,这直接约束了碟形湖对地下水的补给能力。根据先前研究结果,碟形湖以低枯水位季节渗漏补给地下水为主,而高洪水位季节两者之间的交换通量较为微弱,其对地下水的最大补给通量约为0.1 m/d[16]。本文将52个碟形湖做同等处理,均将它们作为源汇项参与到模型计算,根据季节水文变化,设置0~0.1 m/d的动态补给通量。考虑到鄱阳湖洪泛区实际情况,涨水期湖水主要通过地表漫滩(称之为地表路径)以及渗漏补给(称之为地下路径)两种途径来影响洪泛区地下水系统。因此,为进一步考虑溢出面对模型的影响,模型选取丰水期和枯水期两幅遥感影像,通过水面积分布来概化不同时期洪泛区淹水对模型的影响,将其作为地表漫滩的水位边界输入模型(图3a,b)。对于研究区的水文地质参数设定,根据现场钻孔资料(图2b),尽量消除水文地质条件空间差异带来的影响,模型将含水层介质主要作为细砂来考虑。本次地下水模型的时间步长设定为1 d,渗透系数和给水度是本次地下水模型的主要调整参数。为保证非稳定流模拟结果的合理性,模型将多年平均的水文数据作为输入条件,首先开展研究区的稳定流模拟,将模型计算所得的稳定流地下水位作为非稳定流模型的初始条件。本文采用纳什效率系数(Ens)、确定性系数(R2)以及均方根误差(RMSE)对数值模型结果进行定量评价。
图3 鄱阳湖洪泛区枯水期(a)和丰水期(b)水面分布(图中浅绿色表示地表淹水)Fig.3 Water surface area during the dry (a) and flooding (b) seasons of Lake Poyang (the light green color in the figure indicates inundation distributions)
2 结果分析
2.1 模型验证与可靠性评价
本研究中,结合地下水位观测资料,采用手动试错法对地下水模型的主要参数进行调整,调整后渗透系数取值为150 m/d,给水度为0.01。图4绘制了基于野外5个地下水位观测井日水位数据的模型验证效果及定量评价结果,其中4个观测井(1#~4#)位于蚌湖和沙湖湿地周边滩地,1个观测井(5#)位于南矶湿地。由图4a~e可以看出,整个模拟期内,地下水位模拟序列与观测水位序列的年内变化趋势基本一致,尤其体现在地下水位相对较高的时间范围内。然而,模拟序列和观测序列在冬季等一些低水位时期存在一定的偏差,水位偏差大约为0.3~0.4 m。这可能是因为野外地下水位观测井受到降水或周边其它水体汇入的影响,给观测结果带来一定的干扰,导致某些时刻地下水位急剧升高或者下降。其中南矶观测井(5#)冬季模拟水位变化趋势与其它观测井存在一定差异,这可能是由于该观测井位于鄱阳湖上游,其地质条件和水位变化情况不同于下游的几个观测井。通过5个观测点的定量评价结果来看,水位拟合的Ens变化范围介于0.71~0.92之间,确定性系数R2介于0.74~0.94之间,均方根误差RMSE基本小于0.4 m,表明数值模拟结果较为理想,很好地再现了地下水位的时序动态特征以及变化幅度(图4f)。考虑到本文资料有限,且2018年的地下水模型是一个连续、完整的自然年模拟,因此认为当前验证效果表明模型有能力反映地下水的响应变化。总的来说,本文所构建的地下水模型能够很好地适应鄱阳湖洪泛区的水文变化特点,模拟结果的可靠性较高,可以用于洪泛区地下水位的季节性变化对外部环境的综合响应。
图4 洪泛区地下水位模拟效果验证及评价Fig.4 Model validation and performance of simulated groundwater levels in the Lake Poyang floodplain
2.2 地下水动力场特征
2.2.1 地下水位时空分布 图5选取了鄱阳湖“枯-涨-丰-退”4个典型水文时期,以此来分析湖泊洪泛区的地下水位的空间演变规律。考虑到多因素影响下地下水流场响应过程,本文每个典型期选择两个模拟时段加以综合分析。在每年的枯水期,因鄱阳湖主湖区水位较低(约8~9 m),加上研究区受到修水和赣江的地下水入流等补给,整个洪泛湿地通常保持着较高的地下水位(>11 m),总体上地下水由洪泛湿地向东北部主湖区方向流动,即洪泛区地下水补给鄱阳湖主湖区(图5a)。在涨水期,由于鄱阳湖主湖区水位快速提高,因湖水对周边地下水的补给强度发生明显变化,导致主湖区邻近的大多数区域地下水位也随之抬升(约13~15 m),但同时可见洪泛区西北部的地下水位相对较低(<12 m),可知湖水已逐渐向洪泛湿地进行水量传输和补给(图5b)。对于鄱阳湖的高洪水位时期,湖水与整个洪泛湿地融为一体,洪泛区同样受到主湖区水位的强烈影响作用,除了可见邻近主湖区的大部分区域具有较高的地下水位,还可以发现主湖区周边的地下水位要明显高于其它区域,即丰水期仍表现为湖泊补给周边地下水系统(图5c)。退水时期,鄱阳湖主湖区水位逐渐降低(<10 m),但大面积洪泛区仍保持一定的高地下水位(约11~12 m),此时地下水由洪泛洲滩湿地迅速向主湖区排泄,即洪泛地下水补给湖泊(图5d)。上述分析可得,鄱阳湖洪泛区的湖水和地下水位具有明显的季节性转化特征,湖水位的波动变化很大程度上决定了主湖区与周边地下水之间的动态补排模式。总结可知,洪泛区地下水补给湖泊主要发生在鄱阳湖的低枯水位季节(例如枯水和退水阶段),而湖泊补给地下水主要发生在中高水位季节(涨水和洪水阶段)。
图5 鄱阳湖洪泛区枯水期(a)、涨水期(b)、丰水期(c)、退水期(d)的地下水位等值线图(每个水文阶段选择2个代表日期的模拟结果)Fig.5 Contour maps of simulated groundwater levels for dry (a), rising (b), flooding (c), and recession (d) water periods of the Lake Poyang floodplain (simulated results of the two typical dates were selected for different water periods)
此外,通过图6a、b可以清晰地发现,南北方向上,洪泛区地下水位总体上呈现出空间“南高北低”的分布格局,南北区域的地下水水位差可以达到2~4 m,这种空间格局主要与洪泛区地形地貌有关。而在东西方向上,地下水位受主湖区湖泊水位变化控制,涨水和丰水期呈现“东高西低”的变化格局,区域地下水位差异可以达到2~3 m,而退水和枯水期则呈现出“西高东低”的分布格局,区域地下水位差异基本小于2 m。通过图6c洪泛区部分观测点的地下水位变化状况来看,整个洪泛区的地下水位年内波动状况几乎一致,且与湖水位变化规律较为相似。总的来说,鄱阳湖洪泛区的地下水位年内变幅基本介于2~5 m之间。另外,地下水位变幅较大的地方主要位于主湖区附近(例如3#和4#点位),而地下水位变化幅度相对较小的地方则出现在远离主湖区的西侧滩地(例如1#、2#和5#点位)。上述分析表明,就鄱阳湖洪泛区湿地而言,地形地貌对整个洪泛区地下水位分布具有主导作用,但湖泊水位动态变化却是一个关键的外部驱动要素,形成了地下水位时空响应的差异性。
图6 鄱阳湖洪泛区南北和东西断面的地下水位分布(a、b)、典型点位的地下水位变化(c)和剖面位置(d)Fig.6 Distribution of groundwater levels along the south-north and west-east transections (a-b), groundwater level variations for typical observation sites (c), and the corresponding transection within the Lake Poyang floodplain (d)
2.2.2 地下水流速场 为进一步揭示研究区的地下水动力场变化情况,本文选取了不同时期的地下水流速(达西流速)分布结果进行比较分析(图7)。总的来说,整个研究区的地下水流速一般小于1 m/d, 部分区域的流速可达2 m/d,但局部地区因地形起伏变化较大,部分时段可达8~9 m/d。空间上,研究区北部的地下水流速要明显大于南部地区的地下水流速,东部主湖区附近的地下水流速要明显大于洪泛区(图7a~d)。上述地下水流速的空间分布特征主要与地形地貌及其所导致的区域地下水力梯度有关,例如研究区北部的狭窄地形条件以及东部主湖水位对洪泛区地下含水层的激励作用。对于鄱阳湖不同水文时期,地下水流速和运动方向仍存在着明显差异(图7a~d)。枯水季节,地下水主要向东北侧的主湖区流动,主湖区附近的地下水流速较大,洪泛区湿地的地下水流速较小(<1 m/d),在靠近西北侧边界附近的地下水流速则十分微弱(图7a)。涨水季节,地下水流向发生明显转变,由主湖区向洪泛洲滩湿地流动,且由于洪泛区中部地形高程较低,在涨水初期,主湖区和修水、赣江的地下水均向洪泛区中部方向流动,整个空间的地下水最大流速可达1 m/d左右(图7b涨水期选择1时段)。但随着主湖区水位的不断上涨以及湖水-地下水的频繁交换,此时主要表现为地下水整体上由主湖区向洪泛区西部流动(图7b涨水期选择2时段)。在丰水时期,地下水空间流速较涨水期有减小趋势,整体上小于1 m/d,但地下水流动方向仍为由主湖区流向洪泛区,表明了主湖区水位波动的影响作用。因丰水期湖水淹没导致整个研究区相对饱和,大部分区域的地下水流速低于0.1 m/d(图7c)。退水时期,由于东侧湖水位迅速下降,地下水流速的主要方向转变为由洪泛区流向主湖区,地下水流速的空间分布较为明显(图7d)。通过典型断面的地下水流线可清晰发现,在鄱阳湖枯水期和退水期,由于主湖区水位总体上低于洪泛区地下水位,地下水主要由洪泛区流向主湖区,而涨水期和丰水期,地下水由主湖区流向洪泛区(图8)。尽管如此,在研究区南部的部分区域,因地下水位常年低于周边主湖区水位(康山站附近),地下水基本以流向洪泛区为主。上述结果表明,鄱阳湖洪泛区地下水流速场在地形地貌和湖泊水位波动的叠加作用下,流速和流向均呈现明显的季节性差异,这种差异主要体现在地下水在主湖区和洪泛区之间交互作用方式的动态转变。
图7 鄱阳湖洪泛区枯水期(a)、涨水期(b)、丰水期(c)、退水期(d)的地下水流速变化(每个水文阶段选择2个代表日期的模拟结果)Fig.7 Contour maps of simulated groundwater flow velocities for dry (a), rising (b), flooding (c), and recession (d) water periods of the Lake Poyang floodplain (simulated results of the two typical dates were selected for different water periods)
图8 鄱阳湖洪泛区“枯-涨-丰-退”时期典型断面的地下水流线变化(绿线所示的典型断面处为流线起点)Fig.8 Changes of the groundwater streamline for selected transects in the Lake Poyang floodplain during the dry, rising, flooding and recession periods (green lines indicate the start points of the groundwater streamline)
2.3 地下水均衡分析
以整个洪泛区地下水系统为对象,本研究的收入项包括:东侧及地表漫滩一类边界水量补给、降水输入、碟形湖补给地下水量以及西侧二类通量边界补给。主要支出项包括:东侧及地表漫滩一类边界的排泄水量、蒸发量以及西侧二类边界的排泄水量。因此,根据地下水入流量和地下水出流量来计算该地区地下水的储量变化,构建基于月尺度的地下水均衡方程为:
ΔS=DBin-DBout+NBin-NBout+Qp+Ql-Qe
(2)
式中,ΔS为地下水储量的变化量(m3/月);DBin和DBout分别为一类边界的流入量和流出量(m3/月);NBin和NBout分别为二类边界的流入量和流出量(m3/月);Qp为降雨入渗补给量(m3/月);Ql为碟形湖补给地下水量(m3/月);Qe为潜水蒸发量(m3/月)。
根据地下水数值模拟结果,鄱阳湖洪泛区地下水均衡信息汇总如表1所示。由表1可得,2018年研究区的地下水总输入水量约为1.91×108m3,地下水的总输出水量约为1.95×108m3,水均衡计算误差约2%,水量基本保持平衡。根据表1结果分析,降水和蒸发是研究区重要的水量平衡组分,年降水补给量为99.5×106m3,占总补给量的50%以上,地下水的年蒸发排泄量为140.5×106m3,约占总排泄量的70%以上。碟形湖是洪泛区地下水系统的一个稳定补给水源组分,但碟形湖对地下水的年补给量为5.4×106m3,约占总补给量的3%,表明碟形湖对下伏地下水的贡献比重相对较小。另外,第一类边界和第二类边界的年地下水补给量分别为74.8×106和11.2×106m3,占总补给的比重分别为40%和6%左右,同时可见第一类边界的年地下水排泄量(24%)也要大于第二类边界的地下水流出量(4%),表明东侧主湖区对地下水系统(主湖-洪泛区)的贡献作用要明显强于西侧的地下水交换(未控区-洪泛区)。总的来说,鄱阳湖主湖区对研究区地下水平衡的影响要强于修水、赣江等平原地下水等输入作用。
从表1月尺度水均衡变化结果可见,本研究区主要在春、夏季(3-7月)接受外部输入和补给,例如3-5月明显的降雨输入(15.3×106~19.1×106m3/月)以及7月主湖区水量的贡献作用(16.4×106m3/月),而其它月份不同输入组分的贡献相对较小。地下水的蒸发排泄主要发生在夏、秋季(6-10月),月蒸发量变化幅度约为13.6×106~20.4×106m3/月之间。地下水向主湖区的排泄则以秋、冬季为主(特别是12-1月),排泄量总体上大于6.5×106m3/月。图9清晰地呈现了洪泛区地下水系统的主要输入-输出条件和月平均水量动态变化。总体而言,降雨、蒸发和主湖区水量交换是研究区地下水均衡月尺度动态变化的主要组分(图9a),研究区地下水系统在春、夏季以补给状态为主,秋、冬季以排泄状态为主(图9b)。
表1 鄱阳湖洪泛区地下水均衡分析(单位:106 m3/月)*Tab.1 Water budget analysis of the Lake Poyang floodplain groundwater system (unit: 106 m3/month)
图9 鄱阳湖洪泛区地下水均衡组分相对贡献(a)和月平均地下水补(正值)排(负值)量变化(b)Fig.9 Relative contribution of water budget components to the groundwater system (a) and monthly average recharge (positive values) and discharge (negative values) volume (b) for the Lake Poyang floodplain
3 讨论
同其它水源类型相比,地下水是一种极其珍贵的隐藏资源,其水储量要远大于地表淡水河流和湖泊,但地下水的贡献作用极易被忽视[25]。洪泛湿地是生物多样性和生态系统研究的重要区域,具有拦蓄洪水、改善水质等多方面作用,但洪泛湿地的结构和功能对水文变化较为敏感,了解洪泛湿地的地表-地下水转化和其水平衡信息对于探明湿地水文过程、维持生态平衡至关重要[26]。深入认识地下水的组成情况、水量来源及变化特征是研究洪泛湿地可持续发展以及水资源管理的基础保障,更是水质分析、水污染溯源的重要前提[27]。关于洪泛区地下水-地表水交互作用及其重要意义,国外学者已针对亚马逊Amazon湖泊洪泛区[28]、苏格兰Eddleston河流洪泛湿地[29]、澳大利亚Pike和Katarapko河流洪泛区[30]和Murray河流洪泛湿地[31]等开展了大量探索性工作。因此,洪泛湿地的地下水及其动力学机制研究是一个具有前瞻性和挑战性的科学问题。
鄱阳湖复杂的水系统特征、多脉冲的外部输入以及高度动态的水位生消过程,除了湖泊自身地表水动力过程研究,洪泛区地表-地下水转化研究一直是目前众多学者关注的重点。地表-地下水相互作用关系不仅仅是当前水文学、水文地球化学和生态水文学的研究热点,更对其水资源管理、水资源可持续开发利用有着重要意义[3,11]。本文通过分析复杂河湖相互作用下地下水年内动态特征,有效验证了鄱阳湖地下水位与湖水位之间存在转化机制。研究发现,鄱阳湖洪泛湿地的地下水位表现出明显的季节性变化,模拟结果显示,地下水位年内变幅小于5 m,而地表水年内变幅可达9 m以上。由此可见,研究区地下水年内变幅小于地表水,地下水文过程的变化要相对稳定。此外,湖水(低水位出现在第50~60天,高水位出现在第200~210天)与地下水(低水位出现在第60~65天,高水位出现在第205~215天)的高低水位出现时间比较接近,滞后期一般小于15天,说明研究区湖水-地下水文变化的同步性较好。通过进一步分析研究区地下水与湖水位的相互关系可知,越靠近主湖水体的区域,湖水位与地下水位的相关系数越大(可达0.95),而远离主湖水体的洪泛区,地下水与湖水位的相关系数较小(小于0.5)。相关性分析进一步表明,研究区南部的地下水位主要受康山段湖水位影响,而北部地区则主要受都昌段湖水位影响。另外,本文基于单因子分析并采用增加或减小10%、20%偏移量的方法分别对渗透系数K和给水度μ开展了敏感性分析。分析得知,给水度和渗透系数均可以引起较大的地下水位变化,说明模型对给水度和渗透系数的改变都存在着很敏感的响应(图10)。
图10 给水度敏感性分析(a)和渗透系数敏感性分析(b)Fig.10 Sensitivity analysis of the specific yield (a) and the hydraulic conductivity (b)
本文通过FEFLOW数值模拟分析了鄱阳湖洪泛区地下水动力场的演变及地下水的转化通量,提升了对洪泛区地下水系统的理解,但在地下水模型概化过程中仍存在一定的不足:(1)模型的边界条件是影响模型模拟精度的重要因素之一。模型边界主要分为自然边界和人为边界两类。研究区东侧边界为鄱阳湖主河道,可作为自然边界考虑,因此将湖泊按照水位边界进行概化,而研究区西侧很难找到一个完整的自然边界,且研究区西侧主要接受修水、赣江流域地下水补给影响,但地下水流量较小,边界位置对模型的模拟以及水均衡的影响相对有限。此外,本研究区西侧边界为鄱阳湖最大淹没边界,只有在最高湖水位才可能淹没至此,大多情况下该边界对模拟结果产生的影响有限;(2)碟形湖群是鄱阳湖湿地内部的特色水体,模型在考虑碟形湖群对地下水的影响时,将碟形湖的渗漏量作为源汇项进行处理,尽管一定程度上体现了碟形湖群对地下水的影响,但后续工作应充分耦合地下水-碟形湖群的交互作用;(3)地下水数值模型十分依赖于研究区参数给定,同大尺度流域相比,本研究尺度相对有限且参数分区未做考虑。因资料获取有限,模型将洪泛区含水层概化为细砂来考虑,但本文模型最终调整的渗透系数(150 m/d)和给水度(0.01)会与实际背景条件存在一定偏大或偏小的情况,可能很大程度上表征一种砂含砾石的介质特性。实际上,根据水文地质条件的基本特征,洪泛区的非均质性是必然存在的,本文空间均质假定会给地下水均衡组分的计算带来一定的误差。尽管本次模型验证取得了理想的模拟效果,但下一步研究应在大量原位观测和试验基础上,细化研究区含水层结构的空间异质性与水文地质参数分区。同时,针对鄱阳湖洪泛区当前所面临的外部环境变化,下一步工作将围绕湖泊-地下水转化与生态水文的响应,应用气候模式的未来预测结果,考虑鄱阳湖水位极端变化,结合鄱阳湖水利枢纽调度方案,模拟评估不同方案下洪泛区地下水动力场的演变及其对湿地的潜在影响,旨在为湿地生态系统健康的维护、湿地碳中和及未来可持续发展提供支撑。
4 结论
洪泛区地下水对地表河湖系统水资源分配以及生态环境功能等方面具有重要支撑和贡献作用。鄱阳湖洪泛区湿地在长江中下游地区具有重要区位优势和研究特色,但当前变化环境已导致该系统水资源转化和生态演化等朝着不确定性方向发展。本文构建了适应洪泛区特点的FEFLOW地下水流数值模拟,重点开展了洪泛区地表-地下水相互作用及动力学过程研究,综合考虑了洪泛区地形地貌以及气象水文等复杂条件,构建了地下水二维数值模型,模型验证的纳什效率系数变化范围为0.71~0.92,确定性系数为0.74~0.94,均方根误差为0.24~0.43 m,表明该模型能很好地模拟鄱阳湖洪泛湿地的地下水动力特征及对外部环境压力的动态响应,并主要得出以下结论:
1)研究区湖水位的波动变化很大程度上决定了主湖区与周边地下水之间的动态补排模式,即鄱阳湖洪泛区地下水补给湖泊主要发生在枯水和退水时期,地下水主要由洪泛区流向主湖区。而湖泊补给地下水主要发生在涨水和高洪水位时期,地下水由主湖区总体上流向洪泛区。但在研究区最南部,因地下水位常年低于周边主湖区水位,地下水基本以流向洪泛区为主。
2)整个研究区的地下水流速一般小于1~2 m/d,北部地下水流速要明显大于南部地区,东部主湖区附近的地下水流速要明显大于洪泛区。研究区地下水流速场在地形地貌和湖泊水位波动的叠加作用下,呈明显的季节性差异。枯水期、退水期和涨水期的地下水流速较大,基本介于1~2 m/d之间,而丰水期地下水流速较小,大部分区域低于0.1 m/d。
3)地下水均衡分析得出,降雨、蒸发和主湖区水量交换是洪泛区地下水均衡动态变化的主要组分,年降水占总补给量的50%以上,地下水年蒸发排泄量约占总排泄量的70%以上,湖水补给地下水约占总补给量的40%,而地下水补给湖水约占总排泄量的24%。此外,鄱阳湖洪泛区地下水系统主要在春、夏季接受外部降雨输入和主湖区补给,而地下水的蒸发和向湖排泄主要发生在秋、冬季。
4)鄱阳湖主湖区对研究区地下水平衡的影响要强于修水、赣江等平原地下水的输入作用。地形地貌对洪泛区地下水位分布以及流速场具有主导作用,但湖泊水位动态变化却是一个关键的外部驱动要素,形成了地下水位时空响应的差异性以及地下水在主湖区和洪泛区之间交互方式的动态转变。
致谢:感谢南京大学汪祯宸博士,河海大学刘波老师、王哲博士、殷晓然硕士、姜文瑜硕士在模型构建中给予的帮助。