基于数值模型的橡胶坝地下水回补作用研究

李 刚，付信伟

(山东省临沂市兰陵县水利局，山东 兰陵 277700)

0 引言

地下水与河水的相互补给是自然界水循环的重要环节。人造水源工程在防汛供水、蓄水发电的同时，在地下水-河水循环的环节中具有不容忽视的参与度和影响力，对地下水文和生态环境具有深刻的影响[1-2]。人与水的和谐可持续发展是我国现代水利工程建设的核心精神，意为在满足发展需求的限度下最大程度地减少对生态环境的破坏[3]。为此需要对人造水源工程的生态影响效果进行评估。作为一种新型的高分子合成材料的人造水坝，橡胶坝的坝高调节更为灵活，理论上对地下水具有较大的回补作用，但具体的影响数值尚需进一步的研究和检验[4]。此次研究以兰陵县吴坦河橡胶坝水源工程为例，对该区地下水的水位变动数值进行模拟计算，以考察橡胶坝的建设对地下水环境的具体影响。

1 基于数值模型的地下水流场预测

1.1 工程概况及监测点设置

吴坦河属武河行洪道支流的中运河水系。费县和兰陵县交界处的凤凰山是其滥觞[5]。现河道全长44.8km，流域面积达510.97km2，其东西支流汇集于吴坦村，由南至北纵贯兰陵县城，奔向邳苍分洪道。吴坦河地处鲁南郯苍平原，地势北高南低，北部上游为山丘区，出露岩层为寒武系张夏灰岩，以徐庄、毛庄组砂层地层，上部为第四系冲积物所覆盖，土层厚2～10m。中游为河砂，岩层为奥陶系Q1.Q2地层，上部为第四系冲积物覆盖。下游为砂礓盘掺杂，黑褐色黏土覆盖，粘土层厚3m。

北外环橡胶坝工程坝址，位于淮河流域中运河水系吴坦河中游兰陵县城区段，吴坦河中泓桩号26+200处。上游260m处为G206国道吴坦河大桥。工程回水长度4.28km，回水至上游小岭闸。坝址以上流域面积386.43km2，其中包括小马庄水库控制流域面积90km2，长新桥水库控制流域面积36km2。

北外环橡胶坝设于主河槽段，坝高3.50m，净宽70m，共1节坝袋。由坝室段及上下游连接段组成。上游连接段主要有铺盖、岸墙、30m长护坡组成，铺盖为顺水流向15m长，500mm厚C30钢筋砼结构，两侧为C30钢筋砼岸墙；坝室段为顺水流向14.0m长C30钢筋砼底板，厚800mm；下游连接有顺流向长22.8mC30钢筋砼消力池，其中包括7.8m长斜坡段及15m长水平段，消力池深1.4m；消力池下游接顺水流向长15m C25块石砼海漫，两侧为C30钢筋砼岸墙和混凝土联锁块护坡50m；左岸滩地设有泵站，泵站与坝室段经设于坝前铺盖的充排水管路系统相连通。

为监测该区地下水水位的变动，在坝区各坝断面布设6个地下水监测孔，分为6个区域观测吴坦河水位变化、大气降水、其他源汇因素等带来的地下水位变动，以此作为数值模型的计算数据基础。观测孔的孔径为60mm左右，配设自动检测仪，每半小时自动记录观测数据。

1.2 数值模型构建

坝区浅层地下水主要的补给来源是降水和地表水的渗漏，深层地下水主要由上层径流补给。自水源工程建立，地表水的渗漏发生改变，地下水的流场和水位不可避免地受其影响[6-7]。为研究工程对该区地下水的具体影响，可以借助构建数值模型预测地下水流场的状态。这种构建数值模型模拟地下水数值、而后进行计算求解的方法即数值模拟法，是较为通用的地下水模拟方法之一。它有很好的流场分布模拟效果，可以更直观地呈现地下水的水位变化[8-9]。通过模拟地下水数值，可以对地下水流的水量、水位和水质变化进行研究和分析，为评价和管理地下水资源提供帮助[10]。

在假定渗透系数恒定的条件下，采用有限单元法计算渗流数值，首先计算泛函(H)，计算过程如式(1)所示[11]。

(1)

式中，H—总水头，m；σ—渗流的实际阈值；kij—以张量形式表示的渗透系数，m/d；i和j—系数矩阵中的行和列；m—流量增量；s1和s2—边界面。然后如式(2)所示，通过求解等效结点虚流量法计算渗流场的水头分布和河道的渗流量。

[I]{H}={Q0}+{ΔQ1}-{Q2}

(2)

式中，[I]—计算域中总的传导矩阵；{H}—未知的结点水头矩阵；{Q0}—总等效结点的流量矩阵；{Q2}—结点水头已知的条件下，在虚域中对未知水头的等效流量矩阵；{ΔQ1}—结点虚流量矩阵，作用于渗流虚区虚单元。

然后使用二维数值模型模拟该区地下水流的时空关系，模型的表达式如式(3)所示。

(3)

式中，k—渗透系数；(x，y)—坐标变量，m；H—地下水位；H0—初始水头，m；B—含水层地板的标高，m；λ—给水度；i—时间变量，d；D—模拟区域范围；T1、T2和T3—一类边界、二类边界和三类边界；n—二类边界的外法方向；q—二类边界的单位宽度渗流量，m2/d；m′—河床堆积层厚度，m；e—源汇项；Q—河流渗漏补给量，m3/d m2；r—抽水井(注水井)的径向距离，m。模型求解使用开放性和适用性兼备的Modflow软件。

1.3 模型主要参数设置

模型应用之前需设置初始参数值，所需设置的水文地质参数包括基质的渗透系数、给水度、含水率孔隙比等。该区主要岩土层各参数的数值见表1。土层主要成分为壤土，混杂少量中粗砂，分布在河床两岸。壤土的渗透系数为1.5E-05，并且含水率较高，中粗砂的渗透系数为6.8E-03。岩层以溶蚀风化白云质灰岩和裂隙性溶蚀风化白云质灰岩为主，渗透系数分别为1.2E-8～4.6E-9和1.2E-7～7E-10。

表1 各土层主要参数值

设置边界条件时，如式(4)所示计算大气降水的渗入补给量。

(4)

式中，ki—各区域大气降水的渗入补给系数，m3/d；

Ji—各区域大气降水量。m/d；Si—区域面积，m2。

根据实测坝上河道横断面资料量算，坝上水位-蓄水量-水面面积关系见表2。从表中可见，水库蓄水位、库容和水面面积三者正向相关，为模型模拟水库的规模和效益和模拟计算水库渗流提供依据。

表2 北外橡胶坝坝上水位-库容-水面面积关系表

然后使用Modflow的River模块模拟并计算吴坦河地表河水与地下水之间的交换量。运用达西定理分别计算地下水各段的径流补给水量，该部分数据以注水井注水量或抽水井的抽水量的形式引入数值模型。对于潜水的蒸发量，使用软件的Evaportranspiration模块进行模拟计算。该区地下水蒸发的极限深度预估为6m，当超过此界限时，地下水的潜水蒸发量极小，因此忽略补给。蒸发度和蒸发量依据阿维扬诺夫潜水蒸发公式进行计算[12]。对地下水水源的人为开采量数据由历年统计得出，并作为抽水井的抽水量引入所构建的模型。模型参数设定完成后，使用WHS迭代法进行求解。

2 数值模型应用结果分析

模型求解时，WHS算法的最大外部迭代次数为500，最大内部迭代次数为300，算法的收敛精度为0.001m，残差收敛精度为0.001m。阻尼系数是1.0。固定各项渗透系数、边界条件、给水度等参数条件，将模拟的降水、蒸发量、水库水位、潜水开采量等数据输入所建模型，并以2021年10月至2022年9月的水位作为均衡期对照，进行模型验证。验证后的模型利用所用软件的Zone Budget模块均衡分析坝区河水和地下水之间的补给关系，研究橡胶坝的蓄水对该区地下渗流场的作用。模型模拟分析了4种情景下的河水-地下水水位关系。情景一是建坝之后加入河道底部难透水层的情景；情景二是建坝之后不加入河道底部难透水层的情景；情景三是未建坝时加入河道底部难透水层的情景；情景四是未建坝时不加入河道底部难透水层的情景。难透水层虽然渗透性较差，但通过越流也可交换较大的水量，因此作为重要变量纳入考量。

2021年10月模拟水位与实测水位拟合如图1所示，模型模拟的地下水位与实测的水位在南部和西北部部分地区的误差较为明显，模型整体拟合误差较小，可用于进一步的水位预测分析。

图1 模型拟合结果图

根据对地下水水位的模拟计算，得到的四种情景下地下水的均衡分析结果。2021年10月—2022年9月，在情景一下的结果见表3，地下水的总补给量为1582.718万m3，总排泄量为1573.812万m3，均衡差值为8.906万m3，总变化量为8.048万m3。这说明该区地下水流为正均衡状态，满足均衡条件。

表3 2021年10月—2022年9月情景一下地下水模拟均衡分析结果 单位：万m3

在情景二下，地下水模拟均衡分析结果见表4，总地下水得到的总补给量为7856.427万m3，总排泄量为7852.759万m3，均衡差值为3.668万m3，总变化量为2.860万m3。可见，在建坝之后不加入河道底部难透水层时，地下水一定程度上受到了难透水层的影响。

表4 2021年10月—2022年9月情景二下地下水模拟均衡分析结果 单位：万m3

在情景三下，地下水模拟均衡分析结果见表5，总地下水得到的总补给量为1026.821万m3，总排泄量为1020万m3，均衡差值为6.821万m3，总变化量为6.140万m3。从此可知，未建坝时，在模型中加入河道底部难透水层时地下水状态也处于均衡状态。

表5 2021年10月—2022年9月情景三下地下水模拟均衡分析结果 单位：万m3

在情景四下，地下水模拟均衡分析结果见表6，总地下水得到的总补给量为6896.796万m3，总排泄量为6894.451万m3，均衡差值为2.345万m3，总变化量为2.120万m3。未建坝时若忽略河道底部难透水层，坝区河水与地下的水量补给强度较高，二者发生的补给关系较为密切。

表6 2021年10月—2022年9月情景四下地下水模拟均衡分析结果 单位：万m3

总的来看，无论在何种情景下，河流补给是该区地下水补给的主要方式，为其贡献了约过半的补给量。在考虑难透水层的情景下，大坝的建设为该区地下水的水位变动有较大影响。潜水蒸发量与地下水的埋深密切相关，地下水位升高时，蒸发量也会随之增加。在忽略难透水层时，大坝的建设导致了蒸发的加剧。

在4种情景下，坝区河水与该区地下水的水量交换如图2所示。由图可见，除了未建坝时加入河道底部难透水层的情景下秋冬部分月份交换量为负，即发生河流水位较低导致地下水补给河流水的情况，其他情景下交换量均为正。说明在建坝之后无论是否考虑难透水层，河流均对地下水有补给作用，未建坝时不考虑难透水层时河水也对地下水有补给作用。并且，建坝蓄水对该区地下水有较大的水量补给。在情景三下，模拟期的河流为地下水累计补给了36.02万m3的水量。在情景一下，水坝库区预计为地下水累计补给水量为177.21万m3，增加量可达141.19万m3，增幅预计为392.00%。可见，橡胶坝的建成对增加该区地下水水位、回补地下水资源和维护地下水生态具有积极作用。

图2 各情景下河水-地下水逐月交换量

3 结语

橡胶坝工程的建设在改变地表河流面貌的同时，可能对地下水生态造成影响问题。为研究工程对地下水水位变动的影响，研究采用数值模拟法构建了数值模型，模拟了兰陵县吴坦河橡胶坝水源工程对地下水水位变化的影响。模拟分析的结果显示，在诸情景下地下水均处于均衡状态，河流回补了地下水。在考虑河道难透水层的条件下，坝区河流对地下水的预计补给量约为177.21万m3。可见，该模型拟合和预测性能较好，可以为库区地下水水资源开发和水生态保护提供数值分析和理论指导，有利于后续库区水生态评估工作的展开。但研究缺少更充足的观测数据，还需要增加数据来源以提高模型的可靠性。