基于GMS的地下水溶质运移规律研究

张 毅

（乌鲁木齐市水利勘测设计院（有限责任公司），新疆 乌鲁木齐 830049）

0 引言

地下水的流动和溶质在地下含水层中的扩散是一个复杂且不能直接观察的过程，因此如果要具体研究人为抽注水或污染物的排放对含水层的影响就必须建立模型来模拟，从而进行直观的研究。地下水污染问题是环境污染中较严重，而在国内还未引起人们足够重视的问题[1]。地表层下部的构造较为复杂，地下水在地层内的流动比较缓慢，因此，地下水被污染的过程也是缓慢的，其最大的特点就是难以发现、难以治理。一旦地下水被污染，其恢复过程需要数十年[2]。

研究地下水溶质运移规律，就是为了了解污染物在进入地下后的走向，利用数学模型的方法，找地下水及污染物质迁移规律、污染范围以及浓度分布，制定科学的方式来治理地下水污染。本文选定一研究区对其抽水水流运移过程和污染物运移过程进行模拟。

1 研究区背景

本次研究区地下水流数值模拟的范围选定在新疆塔城地区，研究区地形四周高，中间低，整体起伏小，属于中温带干旱和半干旱气候区。

该区含水层厚度约30 m，根据钻孔资料垂直方向上分两层含水层，其中表层岩性为第四系松散沉积物，第二层为砂岩，但在水平方向上水文地质参数不同，可分为三个区域，偏北与偏南区域含水层岩性透水性较弱，中间区域含水层透水性较强。在研究区中央有一抽水井稳定抽水，并在其右侧存在一点源污染及一面源污染[3]。

2 数学模型及数值方法

2.1 数学模型

（1）地下水流模型

用地下水流连续性方程及定结条件来描述地下水流的数学模型。其中初始条件即初始水头由研究区内已知水头的点内插而来，边界条件即研究区北部隔水边界与其余部位定水头边界，并存在一抽水井[4]。

综上，地下水流的数学模型具体如下：

L1，L2分别为研究区北侧及除北侧其余方向的边界。

（2）溶质运移模型

该区存在一点源污染和一面源污染，不考虑化学反应的影响。在地下水流模型的基础上，确定溶质运移的数学模型如下：

I1，I2分别为点源污染与面源污染对应的来源项；D1为面源污染区域。

2.2 数值方法

有限差分法是地下水溶质运移模拟的基础方法。有限差分法就是微分和积分方程的数学解。其基本思想是用有限离散点网格（称为网格节点）代替连续定解区域，用网格上定义的离散变量函数逼近连续定解区域上的连续变量函数；为了逼近原方程和定解条件中的微分商，并用积分和逼近积分，所以原微分方程和定解条件近似地被代数方程组所代替，即有限差分方程组，从而得到离散点上原问题的近似解。然后，利用插值方法，由离散解得到整个区域的定解[5]。

2.3 软件及简介

本次模拟是基于GMS软件进行的，其基本介绍如下[6]：地下水模拟系统（Groundwater Modeling System），简称 GMS，是在已有地下水模型的基础上开发的一个综合性的、用于地下水模拟的图形界面软件[7]。其中MODFLOW为地下水流模拟模块，MT3DMS为溶质运移模拟模块[8]。

3 模型条件

3.1 网格剖分

使用矩形网格剖分，将研究区剖分为一系列立方体，建立网格时，X与Y方向剖分为100×80个，Z方向上设置为2层，见图1。

图1 网格剖分

模拟计算总时长为360 d，时间步长Δt取3 d，定义长度单位为m，时间单位为d，浓度单位为mg/L。

3.2 分层高程

（1）地表高程

研究区地表高程由已知点的高程值内差得，已知点地表高程数据见表1。

表1 研究区已知点地表高程值

续表1

（2）含水层底板高程

研究区分为上下两层含水层，上含水层与下含水层底板高程由已知点的高程值内插而得。

3.3 初始条件及边界条件

3.3.1 地下水流模拟

（1）初始条件：根据给出的已知初始水头点的水头值内差得已知点初始水头值，见表2。

表2 初始水头值

（2）边界条件

研究区北侧为一隔水边界，设置为零流量边界。

其余方向为定水头边界，其中西侧节点定水头为2 m，东侧定水头为1 m，定水头边界根据给出两节点的水头内差而来。上下两层含水层边界条件均为此。

3.3.2 溶质运移模拟

边界条件同地下水流，浓度初始条件为0。

3.4 参数确定

已知研究区地层参数，包括各方向渗透系数（m/d），储水系数，孔隙度以及弥散系数(m2/d)，各分区参数分布见表3。

表3 各分区水文地质参数

3.5 源汇项

1）地下水流模拟

研究区降雨稀少，不考虑降雨入渗与蒸发排泄。

在研究区中部，坐标为（1744.68831，2545.39068）处分布一口抽水井，贯穿上下两层含水层，研究周期内稳定抽水量为350m3/d。

2）溶质运移模拟

研究区中部存在一点源污染（1640.0，2611.5）与一面源污染，污染源贯穿上下两层含水层，污染物为NaCl，不与介质发生化学反应，稳定浓度为5000 mg/L。

3.6 模型识别检验与参数敏感性分析

在抽水井处设置观测点，输入观测水位为0.95，衡量计算值与观测值误差的极差为0.5 m，衡量模型的置信度为95%，将计算结果与实测结果进行对比，可绘出等值线、水头（浓度）-时间变化曲线。见图2。

图2 观测井设置界面

可见散点落在斜率为45°线附近，拟合程度较好。

图3 计算值与观测值拟合图

4 模拟预测结果分析

从研究区大范围看，随着抽水井工作，该区第一含水层地下水流由四周向抽水井汇聚，其中抽水井北侧存在一高水头向四周汇流趋向于水头降低，高水头向远离抽水井方向偏移。而抽水井西侧受抽水影响水流向其汇聚，在保证水源充足的情况下呈上升趋势。总体上，地下水流方向趋向向抽水井汇聚。第二含水层地下水流运动趋势同第一含水层，但第二含水层整体上导水性较强，其响应快于第一含水层。

从抽水井附近水头变化看，第一含水层在抽水井附近的水头变化率较第二含水层大，随着距离的增加，上下含水层的变化趋势接近且水头相近，考虑为第一含水层导水性较差，不利于水流的补给，使抽水井附近的水头变化大，第二含水层导水性强而变化较缓。总体上该区渗透系数较小，抽水井的影响范围较小，见图4。

图4 抽水井附近地下水流模型运行结果

从图4可见，研究区污染物运移范围较小，未抽水时，受地下水流运动影响，污染物向中心（低地势处）运移。抽水井工作后，加速地下水向井处流速，污染物向抽水井方向运移速率增大，在研究期内运移范围扩大。

研究总时长内，点源污染运移范围较面源污染大，结合地下水流模拟结果，抽水井西侧的水力梯度大于东侧，在同一区内，水文地质条件相同，西侧水流速度大于东侧，点源污染物运移速度比面源污染大。而第二含水层导水性较第一含水层大，第二含水层的运移速率较大，相同时间内污染范围较大。该区地下水流方向趋向抽水井，污染源的其余方向的运移速率低，且浓度梯度大，主要受分子扩散影响。

5 结论及建议

该区水流方向本身向低地势汇聚，抽水井工作后改变地下水流方向，形成两个主要的汇水区。同时第一含水层导水性强于第二含水层，在抽水过程中水流补给速度快，而使抽水井附近的水力梯度变化相对较小。该区存在两个污染源，受抽水井影响，污染物趋于向抽水井处运移，且抽水井西侧的水流速度大于东侧，使得点源污染运移速率大于面源污染，其中，第二含水层导水性优于第一含水层，第二含水层的污染物运移速率大于第一含水层。研究区总体上导水性较差，污染物的影响范围较小。

该区污染物影响范围较小，可在其运移速率较大方向设置防渗墙，限制污染物的运移范围，保护地下水免受污染。另外，该区含水层导水系数较小，若抽水井流量过大，水源补给不足，易出现含水层疏干，第一含水层应尤为注意。抽水量结合该区水文地质参数与水资源量确定，该区第一含水层不适应进行抽水，抽水井可设置在导水性较好的region2中地势较低的汇水区，并抽取第二含水层，水量来源更为充足。