基于MODFLOW数值模拟法的地下水饮用水源保护区划分

2013-09-06赵红梅许冠东

四川环境 2013年5期

肖杰，钱骏，赵红梅，许冠东

(四川省环境保护科学研究院，成都 610041)

1 前言

安全的饮用水是人类发展和健康的基本需求，是社会文明进步的重要体现，也是构建和谐社会的基础。目前，我国城乡居民饮水安全状况不容乐观，饮用水水源污染已经成为一个较突出的社会问题，关系到我国的长远可持续发展。据环保部2006～2010年开展的全国城镇集中饮用水水源地的水质调查发现，约1/5水源地存在污染物超标现象，近1/3地下水水源地受到污染。由于我国水污染总体态势未能得到有效遏制，水污染程度仍然在加剧，局部地区的水环境状况恶化，大量因水媒疾病呈频发态势，水污染突发事件屡有发生。因此饮用水源的水质安全保障在现阶段显得尤为紧迫。依据《中华人民共和国水污染防治法》第二十条“省级以上人民政府依法划定生活饮用水地表水源保护区”和《中华人民共和国水污染防治法实施细则》相关规定，设立饮用水水源地保护区，是保护饮用水水源地最大可能免受人类活动影响、保证水质安全的重要措施［1］。

水源地保护区是指为防止水源地污染、保护水源地环境质量而划定并要求加以特殊保护的一定面积的水域和陆域，分为地表水饮用水源保护区和地下水饮用水源保护区两类［2］。随着城市经济建设的逐步发展，土地资源愈加宝贵，水源地保护区划分既要考虑自然因素及水质要求，又要考虑社会经济因素，即水源地保护区的划分在满足保护水质要求的前提下，其面积应尽量小［3］。因此，如何科学合理地划分饮用水源保护区具有重要的研究意义。

2 地下水饮用水源地保护区划分方法

由于地表水水源地保护区划分的模式和经验很多取自于地下水源保护区，因此作者重点探讨地下水饮用水源保护区划分方法。

在国外发达国家中，德国和英国是划分地下水饮用水源地保护区最早的国家，美国起步较晚，但其研究发展十分迅速。德国在水源地保护方面积累了大量经验，其地下水源保护区的划分要求是将取水口所在流域全部划定为水源保护区，在水源保护区内划定3个分区，以取水口中心区保护级别最高，向外逐渐降低，实施分级保护［4］。具体作法为:一级保护区距井周边至少10m，在此区内不能有其它建筑物;二级保护区界以地下水50日流程等值线;三级保护区的范围是水源井降落漏斗的影响范围［5］。随后，英国和美国都形成了自己国家所特有的地下水源保护区划分规范。总体而言，欧美国家划分地下水饮用水源保护区均采用分区保护措施，其水源地保护区的标准也基本类似。

我国对地下水饮用水源保护区的研究起步较晚，水源地保护区的研究经历了以下几个阶段:1984年全国人大颁布《中华人民共和国水污染防治法》并于1996年修正，提出加强水源地保护;1985年国家卫生部发布《生活饮用水卫生标准》规定了水源卫生防护带;2007年国家环保总局发布《饮用水水源保护区划分技术规范》，规定集中供水水源地应设立水源地保护区，其中地下水水源地保护区划分根据不同水源地类型明确了划分原则、方法和要求，将一、二级地下水水源保护区的界线定为100d和1000d流程等值线。至此，我国的饮用水源地保护区研究工作快速发展，国内主要城市都陆续开展了饮用水源保护区划分工作。

近年来，随着研究程度的不断深入，我国地下水饮用水源保护区划分方法主要形成了图形简化法、经验公式法、分析法、圆柱法、水文地质条件法、解析模型法和数值模拟法。其中数值模拟法是求解实际地下水水源地问题的主要方法，采用数值模拟法来计算水源地保护区范围，可以客观并较为详细地刻画实际含水层的结构与水文地质条件，适用于各种背景的水源地研究，尤其对于大型地下水源地保护区的划分，能得到比其他方法更为可靠的结果［6］。

数值模拟法的原理是:利用MODFLOW软件模拟研究区的地下水流场，根据模拟的地下水流场，按不同时间水力截获区的范围来划分水源地保护区。在三维流状态下，t时间截获区是指“在 t时间内，从污染点到水源地的轨迹组所包含的体积”。若忽略垂向流速分量，则截获区为污染物从污染点出发，在t运移时间内，到水源地的水平迹线组所包含的面状区域［7］。对于时间标准的选择，参考《饮用水水源地保护区划分技术规范》［8］中综合我国多地区普遍情况的一般均值，确定水源地一、二级保护区的降解时间分别为:100d、1000d。将污染物运移100d和1000d恰好到达水源地中心的捕获区划分为一级保护区、二级保护区范围，二级保护区以外的地下水补给区划分为准保护区范围。

3 实例分析

崇州市位于四川省岷江中上游川西平原西部边缘，地处东经 103°07'～ 103°49'、北纬 30°30'～30°53'之间，隶属于成都市，东距成都市区37km。崇州市城区集中式饮用水源为地下水，共11口取水井，供水能力 2.5万t/d，水源地距离城区较近，西侧靠近西河，受周边区域污染影响较明显。因此，合理划分地下水饮用水源地保护区对防止水源地污染、保障崇州市的供水安全具有重要的现实意义。本文选取崇州市城区饮用水水源地进行研究。

3.1 水文地质条件

崇州市城区饮用水源地位于崇州市东部平原区，主要分布第四系松散堆积物，包括全新统和中更新统地层。地下水主要赋存于渗透性强的砂砾卵石层中，富水性强。区域地下水补给条件较好，主要接受大气降水和上游地下水径流补给，地下水整体流向为北西～南东向。崇州市城区饮用水源地地下水一部分向西河径流排泄，另一部分向下游径流排泄，地下水水化学类型为HCO3-Ca型水。

全新统 (Q4al+pl)以河流堆积砂、砾石孔隙潜水含水层为主，其二元结构明显，上部为厚0.5～3m粉细砂、粘质砂土，下部为砂砾石层，厚度各地不一。该含水层水量较丰富，补给充沛，渗透性良好，为平原区良好的地下水源，富水程度变化同含水层的分布和规模关系密切。地下水埋藏深度，非汛期1～3m，汛期0.5～1.0m，动态受河水或灌渠水控制，连续降水时影响较明显。地下水化学类型以重碳酸钙水为主，矿化度在0.13～0.40g/L之间。

3.2 水源地开发利用现状及规划

崇州市城区饮用水水源地隶属于崇州市自来水公司。在2008年“5.12”汶川地震后，为有效保证城区人口饮用水安全，崇州市自来水公司经过严格选址与周密筹划、设计，新建11口抽水井抽取地下水供应城区生活用水，该水源地设计供水能力为2.5万m3/d，供水人口11万人，目前，1～9号井已投产，开始向城区供水，现状日供水量2万m3/d。

4 水源地保护区划分

4.1 水源地保护区划分原则

(1)崇州市城区饮用水水源地一级、二级保护区的边界分别为污染质在100d和1000d时间内运移到取水井的捕获区范围，部分边界调整应与实际的水文地质边界相结合。

(2)划定保护区范围应在保证水质达到规定标准的前提下，即确保水源地水质不受污染的前提下，划定的保护区范围应尽可能小，并且便于管理［9］。

4.2 水文地质概念模型

水文地质概念模型是把含水层实际的边界性质、内部结构、渗透性质、水力特征和补给排泄等条件进行概化，便于进行数学与物理模拟。

崇州市城区饮用水源地1～11号抽水井位于西河左岸一级阶地上，地势平坦，地层岩性为河流剥蚀堆积冲积成因的砂卵石，地下水类型为第四系松散岩类 (卵砾石)堆积孔隙潜水。因此，本次模拟含水层对象为砂砾石层。同时根据水源地勘查资料显示，地表下30m内为砂砾石层，没有发现粘性土、粉粒土透镜体。在不影响模拟真实性的前提下，把地表下0～30m视作第一潜水含水层 (Q4al+pl)，30～160m视作第二潜水含水层 (Q3fgl)，两个潜水含水层之间存在垂向补给排泄关系。

本次模型以1～11号抽水井为中心，横向长3360m，纵向长3300m，按照100*100网格划分。区域地表高程为530～555m，取模型高程范围为390～560m，按照垂向上渗透性能的差异，将模型分为3层，局部地区 (抽水井附近)进行网格加密。模型西侧边界为西河河流，可概化为河流定水头边界;东侧边界视作与地下水流向平行，可概括为零通量边界;北西和南西侧可概化为流量边界，其中北西侧为补给边界，而南西侧为排泄边界。模拟区域上边界均为潜水面，在该面上存在大气降水入渗、地表水入渗补给及潜水蒸发排泄等垂向水量交换。模拟区域底部为成都平原第二含水层的隔水顶板－岩性为隔水的泥砾层、粘砂砾石层，该处地下水径流滞缓，上部砂砾石潜水含水层与下界水量交换微弱，基本可忽略不计，因此将其概化为潜水含水层隔水底板。

4.3 地下水水流模型

根据水文地质概念模型建立相应的地下水水流模型［10］。

式中:Kxx，Kyy和 Kzz—分别为渗透系数在 x，y和z方向的分量;

h—水头(L);

W—补给强度(m3/d);

Ss—孔隙介质的贮水率(L－1)或给水率(L－1);

t—时间(T);

Qi—表示区域内第i个井的井下涌水量(m3/d);

δ(x－xi，y－yi)—对应于井坐标的狄拉克函数;

H—表示承压水水头;

Γ1—表示第二类边界(隔水边界);

Γ2—表示第一类边界(定水头边界);

Γ3—表示第三类边界(混合边界)。

图1 示踪粒子运移100d模拟图Fig.1 Simulation of tracer particles with 100 days of migration

4.4 数值模拟结果

建立好水文地质概念模型后，先进行初始渗流场的拟合，对初始水位以及各个参数进行校正。具体做法是:首先对模型进行稳定流计算，将计算结果作为模型的初始水位，再进行非稳定流计算。然后利用实际观测水位和非稳定流模型计算水位进行拟合，在拟合过程中对各个参数进行详细校正。以拟合后的初始渗流场为基础，以抽水井为中心设置环形反向示踪粒子，利用Visual Modflow中MODPATH软件计算出任意时刻质点的位置和流线分布，以示踪粒子逆向运移100d和1000d的范围中最外层质点的起点连接起来就是一级、二级保护区范围。详见图1、图2。