地下水流模型的建立——以华北平原为例

摘要：地下水流模型是辅助地下水资源量评价的一种工具，采用的建模方法主要是有限元法、有限差分法和边界元法。基于有限元法的地下水流模型包括地下水概念模型、数学模型、数值模拟模型3 部分。概念模型建立包括含水层结构概化、研究区边界概化及水力特性确定3 个步骤；数学模型和数值模拟模型建立包括数学方程式描述、区域剖分、定解条件处理、源汇项、水文地质参数确定、模型识别与检验及水量均衡分析7 个步骤。以华北平原地下水流模型为例，阐述建模过程。

关键词：地下水流模型 有限元法 模型识别与检验 均衡分析

中图分类号P641.8 文献标志码B 文章编号1673-4637（2010）02-0018-05

模型是某些复杂的物理实体及其内部流程的一种简化或概化，它能复制或再现一个实际系统的状态[1]，建立地下水流模型的过程就是建立和运行能代表实际含水层行为的模型[2]。地下水模型是一种能帮助分析、解决地下水问题的工具，已经构成现代水文地质学形成和发展的重要推动力之一，成为人们揭示水文地质规律和资源评价与管理中必不可少的工具。地下水模型建模方法主要包括有限元法、有限差分法和边界元法[3]，笔者着重讲述基于有限元法的地下水流模型建立过程，并以华北平原地下水流系统为例。

1 研究区区域概况

1.1 自然地理概况

华北平原在地貌上处于太行山山脉以东，北起燕山山脉，南及东南以黄河为界，东临渤海，整体地势平坦，海拔不超过100 m，自北、西、南3 个方向向渤海湾倾斜。行政区划包括北京、天津、河北3 省的全部平原及河南、山东2 省的黄河以北平原，面积约为13.84 万km2。

1.2 水文地质概况

华北平原第四系是由多层砂、砾石和黏土交叠、纵横交错构成。从山前平原、中部平原至滨海平原砂砾石层，颗粒由卵砾石、中粗砂到细砂、粉砂，厚度在山前平原由薄变厚，至中部平原边部变薄，复而沿沉积方向又加厚，至滨海平原又逐渐变薄。

2 地下水流概念模型建立

在分析研究区水文地质条件的基础上，建立研究区三维水文地质概念模型，包括含水层结构概化、边界条件概化和水力特性确定。

2.1 含水层结构概化

结合含水层组的划分，概念模型中，将含水岩组在垂向上概化为3 个含水层。

2.2 边界概化

2.2.1侧向边界

研究区西部和北部山区与平原自然分界线概化为流量边界，模型第1层经过此边界接受山区侧向补给；模型深层的山前边界在本次计算中概化为隔水边界。南部及东南部以黄河为边界，模型第1层与黄河有水力联系，概化为流量边界；其他层位概化为隔水边界。东部以渤海海岸线为边界，模型第1 层概化为水位边界，下伏其他各层定为隔水边界。

2.2.2垂向边界

潜水含水层自由水面为系统的上边界，通过该边界，潜水与系统外界发生垂向水量交换，如接受大气降水入渗补给、蒸发排泄等。第1 含水层和第2 含水层通过微弱的越流交换物质和能量，其越流量由相邻两层的水位差及垂向渗透系数、含水层厚度决定。研究区底部边界概化为隔水边界。

2.3 水力特性确定

从空间上看，研究区地下水流整体上以水平运动为主、垂向运动为辅，地下水系统符合质量守恒定律和能量守恒定律；在常温常压下地下水运动符合达西定律；考虑2 个含水层之间的流量交换，地下水运动可以概化为空间三维流；地下水系统的垂向运动是由层间水头差异引起的；地下水系统的输入输出随时间、空间变化，故地下水为非稳定流；参数随空间变化，体现了系统的非均质性，所以含水介质概化为非均质各向同性介质。

综上所述，研究区可概化成非均质、水平方向各向同性、垂向存在变异、空间三维结构、非稳定地下水流系统，即地下水系统的概念模型。

3 地下水流数学模型和模拟模型

3.1 地下水流数学模型

对于上述非均质、各向同性、空间三维结构、非稳定地下水流系统，可用地下水流连续性方程及其定解条件来描述。选择地下水模型软件FEFLOW 求解该定解问题，以建立研究区地下水数值模拟模型。

3.2 地下水流数值模拟模型

研究区地下水流数值模拟软件采用基于有限元法的FEFLOW 软件，FEFLOW（Finite Element subsurface Flow system）是由德国WASY 公司开发的，它是迄今为止功能最为齐全的地下水模拟软件包之一，可用于复杂三维非稳定水流和污染物运移的模拟[4]。

3.2.1网格剖分

采用不规则三角剖分，剖分后的模拟区共有8 万3 805 个结点，12 万8 948 个单元格，其中每1 层为近2 万7 935 个结点，对于华北平原约14×104 km2 的面积上，约5 km2 为1 个网格。

3.2.2定解条件处理

3.2.2.1 初始条件

根据2001 年12月研究区内统测的地下水位监测资料，按照内插法和外推法获得各层初始水位。

3.2.2.2 边界条件

整个研究区是一个相对完整的地下水系统，西部及北部的山区和平原区之间的自然分界线作为侧向流入边界；南及东南部黄河处理为流量边界；东部渤海边界处理为一类水头边界。

3.2.2.3 模拟期

模拟时期为2002 年1月到2003年12月，为2 个完整的日历年。以1个月作为1个应力期（即模型运行的时段：在相同的时段内，各项应力包括补给项和排泄项等不发生变化）。每个时间段内包括若干时间步长，时间步长为模型自动控制，根据模型运行情况来自动调整时间步长长度，并严格控制每次迭代的误差。

3.3 源汇项处理

地下水系统的源汇项是指其补给项和排泄项。在均衡区内地下水补给项主要包括有降水入渗补给量、农业回灌量、侧向补给量、河流渗漏量；地下水排泄项主要包括开采量和潜水蒸发量。源汇项处理的目的在于确定地下水的各个补排项随时间和空间的变化规律，为建立地下水流数值模拟模型准备数据。

3.3.1补给项

（1）降水入渗补给。

本次计算是在充分考虑不同降水量、不同水位埋深以及包气带岩性等因素的情况下进行的。结合资料特点将降水入渗系数按面积加权平均到每个行政县市。

在模型中再将各小区的降水入渗量处理为各行政区县上的降水入渗强度即单位面积上的降水入渗量。全区总量2002、2003 年2 年分别为110.64×108m3 / a、213.82×108 m3 / a。

（2）农业回灌量。

农田回灌量主要取决于灌溉定额、灌溉次数以及耕地面积和包气带岩性。研究区有3 种灌溉方式：井灌、渠灌和井渠结合。

全区总量2002、2003 年2 年分别为37.27×108m3/a、35.6×108 m3/a。

（3）侧向补给量。

侧向补给量包括华北平原山前侧向量和黄河侧渗量。

计算结果: 2002、2003 年华北平原山前浅层侧向补给量分别为14.26×108 m3 / a、14.79×108 m3 / a，深层侧向补给量均为5.24×108 m3 / a； 2002、2003 年黄河侧渗量分别为3.25×108 m3 / a， 3.19×108 m3 / a。

（4）河流渗漏量。

河流渗漏量是在对河水与地下水的补排关系分析后，确定出河水补给地下水的河段，采用地下水动力学法。计算结果: 2002 年、2003 年2 年华北平原5个省市的河流渗漏总量分别为2.53×108 m3 / a、2.56×108 m3 / a。

3.3.2排泄项

（1）开采量。研究区内开采量统计到各个区县。开采量包括农业开采量、生活和工业开采量以及少量的生态用水。经计算，2002、2003 年华北平原浅层地下水开采量为176.62×108 m3 / a、172.63×108 m3 / a，深层地下水开采量为31.24×108 m3 / a、30.8×108 m3 / a。

（2）潜水蒸发量。潜水蒸发量主要与潜水位的埋深、包气带岩性、地表植被和气候等因素有关。根据前人的研究成果，一般认为该区的蒸发极限埋深为（4～5）m。利用阿维扬诺夫公式来计算蒸发量。经计算，全区2002、2003 年年蒸发总量分别为62.53×108 / a、50.19×108 / a。

3.4 水文地质参数确定

水文地质参数是表征含水介质储水能力、释水能力和地下水运动能力的指标。对于华北平原这个典型的大区域，水文地质条件极其复杂，水文地质参数分区各含水层组差别较大，为了能更详细地刻画研究区水文地质条件，给出了不同含水层的参数值范围，第一含水层主要的水文地质参数为水平渗透系数和包气带的给水度；第二含水层、第三含水层主要的水文地质参数为渗透系数、释水系数。

3.5 模型识别与检验

地下水模型的识别与检验过程是整个模拟中极为重要的一步工作，通常要进行反复地修改参数和调整某些源汇项才能达到较为理想的拟合结果，同时也是对水文地质条件再认识的过程。它直接关系到后续的地下水资源评价、水位预测和科学管理的可靠性与精度问题。3.5.1模型的识别和验证遵循原则

（1）模拟的地下水流场要与实际地下水流场基本一致，即要求地下水模拟等值线与实测地下水位等值线形状相似。

（2）模拟地下水的动态过程要与实测的动态过程基本相似，即要求模拟与实际地下水位过程线形状相似。

（3）从均衡的角度出发，模拟的地下水均衡变化与实际要基本相符。

（4）识别的水文地质参数要符合实际水文地质条件。

3.5.2水位拟合

根据以上4 个原则，通过反复调整参数和均衡量，识别水文地质条件，确定了模型结构、参数和均衡要素。模拟期末（2003 年12 月31 日）浅层、深层地下水的模拟流场与实际流场对比见图1、图2。浅层、深层地下水典型长观孔水位过程线拟合情况见图3。

由图1，图2 可知，模拟潜水流场与实际流场总体流动方向相同，从过程拟合曲线可以看出模拟的地下水动态和实际动态变化趋势基本一致。模拟研究区流场与实际流场总体流动方向相同，主要漏斗区范围基本一致；体现了区内地下水的运动规律即整体趋势从山前向中部平原再到滨海平原，局部人为开采造成的漏斗区形成局部流动系统。识别后的水文地质参数见表1。

表1　   识别后的水文地质参数表

3.6 水量均衡分析

通过对地下水各均衡要素分析，得出模拟区地下水系统浅层和深层含水组在模拟期2002 年1 月至2003 年12 月里水量均衡结果，见表2。

表2　   模拟期水量均衡分析表

研究区在模拟期内地下水系统总补给量为465.66×108m3（不包含层间的越流量），总排泄量522.97×108m3，均衡差为-63.27×108 m3，为负均衡。模拟期2002 和2003年，浅层含水组总补给资源量为449.75×108 m3，降水入渗量为336.30 ×108 m3，占总补给资源量的74.77％，灌溉回渗量为72.87×108 m3，占总补给资源量的16.20％，山前侧向流入量约为29.05×108 m3，占总补给资源量6.46%，黄河侧渗补给量与河道渗漏量分别占1.43％和1.13％。在总补给量中，降水入渗补给量占绝对地位，其次是灌溉回渗和山前侧向补给量。

模拟期内浅层含水组总排泄量为490.41×108 m3，人工开采为349.25×108 m3，占总排泄量的71.22%，蒸发为112.72×108 m3，占总排泄量的22.98％，其次层间越流排泄量为28.44×108 m3，占总排泄量的5.8％。在所有排泄项中，人工开采所占的比重接近3/4，为地下水主要的排泄途径，其次为蒸发排泄。

3.7 小结

综上所述，从模拟期末等值线拟合、典型长观孔水位过程线拟合以及研究区地下水系统均衡分析来看，通过水文地质结构模型结合水文地质概念模型建立的研究区地下水流数值模拟模型基本达到精度要求，符合研究区实际的水文地质条件，也基本上较好地反映了研究区内地下水的运动特点和动态变化特征，故可利用模型进行地下水资源评价和开采条件下地下水流场演化的趋势性预测。

4 结语

（1）地下水模型是一种能帮助分析、解决地下水问题的工具，可用来进行水资源评价及各类开采条件下，地下水流演化趋势预测。

（2）地下水模型的识别与验证较前期建模工作复杂得多，需要反复修改参数和均衡项，以求得模型能准确反映实际水文地质条件。

（3）模型精度视研究区的范围、复杂程度、资料精度以及用途需求而定，不能用水源地模型的精度范围来要求大区域地下水模型的精度。

参考文献

[1] 孙纳正. 地下水流的数学模型和数值方法[M]. 地质出版社，1981.

[2] J ames W . Mercer，Charles R. Faust . Ground -Water Modeling，National Water Well Association[M]. 1981.

[3] 张宏仁，李俊亭. 有限差分法与有限单元法在渗流问题中的对比[J]. 水文地质工程地质，1979（2）.

[4] 孙祥光，王井泉，翁明华. Hans-J觟rg G D. FEFLOW 有限单元地下水流系统[M]. 中国矿业大学出版社，2000.1-154.

作者简介：林坜（1980—），女，工程师。