基于GIS的地下水模型边界地下水流估算方法

于国宝

(辽宁省朝阳水文局,辽宁 朝阳 122000)

1 概述

地下水和河流系统的生态完整性经常受到人类活动的威胁[1]。社会也面临日益严重的水资源管理问题[2]。通常地下水建模是确定地下水管理策略的重要工具，然而，为了开发一个可靠的地下水模型，研究人员需要克服许多方法上的挑战[3]。本研究论文集中于一种基于地理信息系统(GIS)的新方法，用于在地下水流模拟中估算研究区域边界上的地下水流。

作为开发地下水流动模型的第一步，需要确定研究区域和模型边界。模型边界是研究区域和周围环境之间的界面[4]。阿慧娟等[5]证明了分配正确边界条件的重要性，并表明边界位置的不准确性会导致校准过程中引入不必要的不均匀性。同时河流、湖泊和湿地与地下水的相互作用受水体相对于地下水流动系统的位置、河床的地质特征及其气候环境的控制[6]。陈张建[7]全面总结了河流和地下水之间的相互作用。并认为在降雨量少的情况下，许多溪流中的基流大部分时间是来自于地下水的输入。

GIS是为任何地下水流和污染物迁移建模工作开发地下水模型的重要工具。GIS提供数据管理和空间分析能力，可用于地下水建模。它提供数据插值、系统的模型参数分配、空间统计生成和模型结果的可视化显示[8]。唐敏等[9]设计了一个GIS数据库，研究了地下水脆弱性分析和水文地质建模。GIS在地下水流系统的概念化、特征化和数值模拟中也具有广泛应用。庞园等[10]应用GIS技术对坚硬岩石地形的地下水补给进行定量建模。而且GIS还提供了许多工具，如创建剖面图的工具、提取某一点的值、创建缓冲区的工具、插值工具等。因此，它可以和达西定律一起用来估算通过给定边界的平均流量。

2 研究区

本研究选取辽宁省浑河的某支流作为研究对象。研究区的主要含水地层是冲积层，约占总面积的60%，其次是花岗岩和千枚岩。研究区的地下水位深度从平均海平面(MSL)的136～380m不等。MSL的地面高度从147～446m不等。该地区的平均降雨量为452mm。大部分降雨发生在6月至9月的季风季节。其余时间被认为是非季风季节。本文从水文、水文地质、地下水流方向和形态特征等方面对该区域进行了深入研究，并在此基础上确定了地下水流的边界条件。在研究区域的东侧，是倾向NNE-SSW的山脉。可以将其视为一个重要的无流量边界。在西侧剩余的小部分以封闭研究区域的边界。研究区域如图1所示。选择了研究区域内及其周围的108个探井，如图2所示，根据这些探井的历史水位数据，用逆距离加权(IDW)插值法计算了1999—2003年季风和非季风季节的水位面。插值时选择360m×360m的像元尺寸作为输出图像。图2还显示了2003年季风前地下水位的分布图像。

图1 研究区域

图2 测点位置

根据获得的28个探井记录，绘制了研究区的水文地层图。建造了3个冲积层、花岗岩和千枚岩的水文地层，并计算出所有三层的底部高程。ASTER高程值用于开发数字地形模型。研究区部分探井的详细信息以及所有三层的井底高程见表1。利用钻孔测井信息和IDW插值工具，为冲积层、花岗岩和千枚岩的水文地层底部高程准备表层。还计算了监测井所有季节的平均水位深度，并通过插值将其近似转换为1个层。

表1 部分探井数据

3 研究方法与结果讨论

本文采用了以下方法对边界进行计算。

(1)在研究区的边界上创建了一个5公里的缓冲区。边界被分成22个横截面，如图3所示。西侧的三个横截面编号为F1、F2和F3。上侧的部分从U1到U10顺时针编号，底端的部分从L1到L9顺时针编号。图3显示了其中一些部分的编号。在遇到边界急剧弯曲的地方，应采用更多的横截面。

图3 横截面示意图

(2)利用三维分析工具，从1999—2003年的不同季风前和季风后水面剖面中提取了沿每个横截面区域大约30—35个均匀分布点的地下水位。

(3)绘制这些点，并将四阶多项式拟合到每条曲线，以确定精确边界点的梯度，如图4所示。观察到四阶多项式的适合度在所有情况下都非常好。2002年后季风的U4截面拟合相关性为0.998。从四阶方程出发，通过简单的微分，可以更容易地得到任意点的梯度(dy/dx)。因此，在每个这样的轮廓的边界点计算梯度。

图4 2002年季风后沿U4断面的水位变化

(4)尽管不同剖面之间的梯度值不同，但观察到季风季节前后和不同年份之间给定剖面的梯度计算值没有显著变化。典型线段几个横截面不同年份的季风和非季风季节的梯度计算值见表2。流入研究区域边界的梯度值保持为正，流出为负。可以明显的发现，在季风前后各年份和各位置地下水位并没有很明显的变化，特别是对于L9而言，水位梯度变化值小于0.00001。

(5)因此，可以通过取所有季节所有梯度计算

表2 典型线段的梯度计算值

值的平均值来计算每个横截面的梯度平均值。然后对所有横截面获得的边界段的平均梯度值进行分析，并识别和合并具有几乎相同梯度值的相邻边界段。通过这种方法，整个边界被重新分组为8个边界段，并计算梯度的平均值。对于这8个边界段，基于单个线段的平均值。这些分段和梯度平均值如图5所示。东部边界由于靠近山脉，因此显示值为“0”(无流量边界)。

(6)在每个边界段选择5个点，并从相应的层中提取这些点的地面高度、水文地层底部高度和平均水深值。因此，在8个边界段上选择的40个边界点的细节制成表格。然后用这些值计算三层含水层系统在每个点的等效横向渗透系数，然后用达西定律计算各点边界每米长度的净流量，再计算每个区段的平均流量值，并乘以区段长度，计算每个区段的总净流量，见表3，大部分的梯度值很低。梯度最大值(0.0024)发生在在4号段。用达西定律计算得到的流量结果表明单位长度下的净流量比从-0.5040到0.9778之间变化，在第2段流入量最大，为35725m3/s。在第6段流出量最大，为29000m3/s。

表3 每个边界段的净流量

4 结论

本文在无法确定研究区域附近无流动边界的情况下，构建了一种基于地理信息系统计算研究区域边界的净流量的方法。并认为地下水位梯度在不同季节和不同年份之间不会发生显著变化，此假设也在地下水模型中得到验证。在本研究中，研究区域边界分为9个区段，并计算每个边界区段的流量梯度值。通常，梯度值很低。梯度最大值(0.0024)发生在在4号段。然后使用达西定律来计算穿过这些边界的流量。结果表明：在第2段流入量最大，为35725m3/s；在第6段流出量最大，为29000m3/s。本文提供了一种计算边界平均净流量的方法，以期对地下水污染研究和地下水模型建立提供一定的参考。

图5 不同边界段的最终水位梯度