第5章 华北平原非饱和带中地下水滞留时间

包气带中地下水运动时间（地下水自地表运动至潜水面所需的时间）是影响浅层地下水年龄，评价地下水是否受到地表污染影响的重要因素。点上的包气带中的地下水运动时间可通过示踪剂的运移模拟估算。计算示踪剂进入包气带到达潜水位的运动时间作为地下水在包气带中的运动时间。而区域尺度上非饱和带中地下水运动时间可通过地下水补给量和非饱和带的水文地质特性估算。

5.1 区域地下水补给计算方法

在干旱、半干旱地区,农业依赖于灌溉,因此地下水补给以及蒸发量的准确估算对实现地下水资源的可持续管理极为重要。然而,由于土壤质地和植被类型的空间变化、以及气象因素,如降水和蒸发的时空变异,地下水补给量通常不易估算,并且在灌溉的条件下,地下水的补给过程会更为复杂(Gee和Hillel,1988)。传统的测量地下水补给的物理方法(如蒸渗仪和渗透仪lysimeter,infiltrometer)和测量蒸散发的物理方法(如涡度相关和蒸渗仪),以及诸多经验方法通常要求大量的试验数据,因此往往局限在小尺度应用。区域地下水补给量的估算通常采用的方法有:水均衡法、潜水位动态法(WTF)、化学和同位素方法以及数值模拟方法(Simmers,1987;Healy和Cook,2002;Scanlon等,2002;Xu和Beekman,2003)。仅应用某种方法估算地下水补给量的准确度一般很难判定,因此建议同时采用多种方法以对不同方法的结果互相验证以降低结果的不确定性(Scanlon 等,2002)。另外,由于气象因素、土地利用类型以及地表物理和生物特性的空间变异,将点状计算结果插值到区域上往往造成较大误差(Scanlon等,2002)。因此,对地下水补给估算方法的选择取决于研究的空间和时间尺度,以及研究区的数据情况(Scanlon等,2002;Delin等,2007)。 估算区域地下水补给量的主要方法为土壤水量平衡法,其中,土壤带或根植层底界的地下水流出量作为净入渗量(net infiltration)或潜在补给量(potential recharge)。基于该种方法的数值模型很多,例如美国地质调查局(USGS)的INFIL 模型(USGS,2008),BALAN 模型(Samper和García-Vera,1992),WAVES模型(Zhang和Dawes,1998),和美国环保署(U.S. EPA)的HELP模型(Schroeder 等,1994)。虽然这些模型能够简便快速的估算地下水补给量的区域分布,前提条件是地下水潜水面要低于根植层的底界,并假设根植层底部的渗出量即为地下水补给量(潜在补给量=实际补给量),认为自根植层底界至潜水面的非饱和带内的水流过程是可以忽略的。然而,当非饱和带,或包气带厚度较大时,非饱和带将成为控制补给时间和水分传输的主要因素,忽略其中的水流过程是不可取的 (Harter和Hopmans,2004;Hunt等,2008)。 鉴于以上各种方法的优缺点,非饱和带数值模拟方法逐渐成为估计区域地下水补给时空分布的通用方法(例如,Keese等,2005;Small,2005;Liggett等,2009)。此类模型主要包括 DAISY(Hansen等,1990),UNSAT-H(Fayer,2000),SWAP(Kroes和 van Dam,2003)和HYDRUS -1D(Šimůnek等,2005)。虽然这些模型能够模拟非饱和带中的水流并估算蒸发量和地下水补给量,其应用过程中通常需要指定模型下边界条件(如潜水面位置),没有耦合饱和带地下水流动过程。因此其结果不能直接反映出气候变化、土地利用和土壤质地对地下水动态的影响。近年来开发的耦合非饱和带和饱和带水流过程的数值模型,如HYDRUS-3D(Šimůnek,2006),MODFLOW-VSF(Thom等,2006)和HydroGeoSphere(Therrien等,2005),是对上述非饱和带数值模型的巨大改进。然而,这些模型的理论基础是刻画非饱和带水流过程的理查德(Richards’)方程,其应用主要受到计算机计算能力的限制(van Walsum和Groenendijk,2008)。而对于流域尺度的模拟计算,将非饱和带内水分运动假设为一维流动通常是合理的(Harter和Hopmans,2004)。因此,许多模型,如MIKE-SHE(Refsgaard等,1995)和MODFLOW-UZF(Niswonger等,2006)采用了非饱和带一维流动和饱和带三维流动的耦合方式以降低耦合模型的复杂程度和提高计算效率。在这些模型中,一维的理查德(Richards’)方程被进一步近似成运动波(kinematic-wave)方程并采用特征值方法求解。该求解方法的优点在于不需要对非饱和带进行网格剖分(Niswonger等,2006)。 基于非饱和带水流过程模拟估算区域地下水补给量的主要困难在于模型参数的获取和尺度转化,如土壤的渗透系数。对于该过程,Harter和Hopman(s2004) 对多种地质统计方法提升土壤渗透系数的空间尺度做了细致的综述性讨论,但这些方法是基于大量的小尺度范围内的渗透系数估算结果(Vereecken等,2007)。作为直接测量原状土壤渗透系数的替换方法,土壤转换函数方法(PTFs)(Wösten,2001)能有效的用于区域非饱和带渗透系数的估算。利用土壤转换函数方法可以直接测量土壤特征指标(通常附带于土壤调查数据库)并估算非饱和的水力参数(Li等,2007)。许多广泛应用的土壤转换函数方法是基于统计模型的方法(Pachepsky和Rawls,2004)估算土壤水特征函数和饱和渗透系数。所利用的土壤特征指标通常为土壤颗粒分级,土壤干密度和有机成分含量(Nemes等,2003)。

绝大多数非饱和带数值模拟软件最初的开发目的是用于农业领域研究工作,需要大量反应植被生长特征的参数。虽然土地利用和植被类型分布可以通过调查统计直接获得,遥感解译也有助于获取植被参数,但对于植被类型和土地利用不是地下水补给主控因素的地区,简化植被生长过程模拟,减少模型参数更有利于减少模型的不确定性和模型校正过程中的不唯一性。 蒸散发(ET),作为土壤水均衡的主要排泄项,对利用包气带数值模拟方法计算地下水补给有着重要影响(环境同位素方法估计地下水补给不需要蒸散发的计算)。上述不同的模型中对蒸散发的计算过程也各不相同。最常用的计算方式首先是基于常规气象观测数据计算潜在蒸发或参考蒸散发(ETp或ET0)(Lu等,2005),然后利用经验公式计算裸土蒸发和植物蒸腾(如Šimůnek等,2005)。实际蒸发和植物蒸腾根据土壤含水量情况和植被信息确定,并加和得到实际蒸散发(ETa)。但是,这种方法通常要求较小的时间尺度,如日气象数据。因此,饱和带水流模型的时间步长也需要相应减小,但对于饱和带地下水流、特别是大区域地下水流模型来说,显然是不合适的。作为一种折衷方案,可以利用子时间步长,即在饱和带水流模型时间步长的基础上再进行剖分,得到的更小时间步长作为非饱和带模型的时间步长,如耦合HYDRUS和MODFLOW的程序包(Seo等,2007)。 基于互补关系的蒸散发模型(Bouchet,1963)是另外一类用于估算区域蒸散发的常用模型,这类模型仅利用常规气象观测数据,不需要刻画陆面边界的详细信息。因此,此类模型避免了复杂的土壤-植被-大气这一复杂系统,不需要土壤含水量数据和植被参数(Hobbins等,2001)。不同的气候类型地区与区域水均衡的蒸发计算结果对比已经证实了该类方法的适用性(如,Hobbins等,2001;Xu和Singh,2005)。此外,此类模型认为陆表大气湿度包括了所有陆面过程的影响,已经考虑了灌溉对蒸发的影响(例如,Ozdogan和Salvucci,2004)。 

5.2 华北平原地下水补给研究现状

前人研究认为华北平原降水入渗补给占到地下水自然补给量的70-80 % (陈望和,1999;陈梦熊和马凤山,2002)。华北平原地下水资源评价工作中通常采用渗透系数方法估算地下水补给量,入渗系数一般采用均衡试验得到,主要受土壤岩性控制,与降水量无关。计算得到的区域地下水补给量范围较大,为 70-180 mm/yr(Zhu和Zheng,1983;刘昌明和魏忠义,1989;陈望和,1999;段永侯和肖国强,2003;任宪韶,2007)。该方法是对地下水补给过程的高度简化,相当于认为地下水补给和降水量或灌溉量之间是一种线性关系(例如,Kendy 等,2003),然而实际的补给系数受土壤特性、植被、包气带厚度和气象条件等多种因素控制,且所有这些控制因素都是随时间和空间变化的(Scanlon等,2002;Li,2009)。 应用环境同位素方法 (Wang等,2008;Liu等,2009;Song等,2009;von Rohden等,2010;Yuan等,2011)和土壤水均衡和非饱和带数值模拟方法(例如,Kendy等,2004;Lu等,2011)可以估算小尺度地下水补给量。这些小尺度地下水补给量估算结果在量级上基本一致(100s mm/yr;见表 5-2),这些计算结果仅能代表研究场地一定范围内特征地下水补给模式(降水入渗补给或灌溉回归补给),并与所采用的补给方法有关(见表 5-2)。华北平原缺乏对集中式补给的定量研究,但这些补给仅限于局部地区(如,Yuan 等,2011)。鉴于此,本研究进行了潜水面的地下水补给量估算。 影响华北平原地下水补给的一个关键因素是非饱和带厚度的增大。华北平原长期大量开采地下水造成包气带厚度从70年代的2-15 m 增大到8-30 m甚至局部地区达到30-56 m(张光辉等,2004)。忽略非饱和带中的水流过程,以及非饱和带厚度对地下水补给的影响,特别是对地下水污染的评估和地下水保护策略的制定都是不可取的。虽然已有研究利用土壤含水量和土水势观测评估非饱和带厚度对入渗过程的影响(张光辉,1992;张光辉等,2007),但很少涉及非饱和带对区域地下水补给的影响研究。

另一个关键问题是如何将局部地下水补给估算结果和模型参数提升到区域尺度,即水资源方案制定层面必须考虑的尺度范围。一般来说,包括非饱和带厚度和渗透系数,降水模型,灌溉类型及模式,土地利用以及蒸散发在内等大量数据是提升补给量计算尺度的必要条件。然而,考虑到华北平原的空间尺度及现有数据和研究成果,利用局部地下水补给量估算区域地下水补给量显然是不现实的。虽然最近有研究利用数据融合技术耦合模型计算结果和环境同位素信息计算结果来估算区域地下水分布(如,Ng等,2009),但现阶段来说对区域研究计算量仍很大。本次研究中采用的方法对区域地下水补给估算是适用的,并能够考虑控制地下水补给过程的一些关键因素。  本研究耦合土壤水均衡模型与MODFLOW的UZF程序包,建立区域地下水流数值模型,模型能够模拟非饱和带水流过程并估算地下水补给量。该方法能够利用多种有效估算区域陆面蒸散发的方法,如遥感方法,并且适用于区域尺度的模型校正和应用。本研究的主要目的是评价土壤质地类型,非饱和带厚度变化对华北平原区域地下水补给的影响。华北平原的空间尺度,土壤转换函数和已有数据能够支持建立合理适用的非饱和带-饱和带耦合水流模型来模拟地下水的动态演变。模型校正不是本章研究的主要目的,但通过对比模型计算的地下水补给量和与其他估算方法取得的结果,以及计算地下水位和观测地下水位来验证计算得到的地下水补给量(见表 5-2)。 

5.3 华北平原非饱和带-饱和带耦合数值模型

5.3.1 土壤水均衡模型

当没有明显地表径流时,某特定时段内根植层的土壤水均衡可表示为:

式中Rp为根植层底部地下水渗出量,或潜在地下水补给;P为降水;Ir为灌溉; ETa为实际蒸发量; 为土壤水储存量变化量。

5.3.2 非饱和带水流模型

描述均质非饱和带中垂向一维水流过程的理查德方程可表示为: 

式中θ土壤水体积含水量;S为单位长度的蒸散量;K(θ)为渗透系数函数,与含水量有关;z为垂向上高程;t为时间。 MODFLOW的UZF非饱和带程序包(Niswonger 等,2006)对公式(5-2) 进行了简化,去掉了扩散项,假设垂向上的地下水流动仅有重力驱动。UZF中的控制方程表示为:

一般来说,基于理查德方程的非饱和带水流模型需要较细的网格剖分,因此对于计算时段较长的模型来说,要求较高的计算能力(Hunt 等,2008)。再者,建立非饱和带渗透系数与土壤含水量关系函数的参数在区域范围内不能直接获得。对于区域尺度上的非饱和带水流模拟,能够获得控制方程的闭合形式解,并保证数值计算的稳定性,因此更有利于大空间和长时间尺度的模型建立。鉴于此,UZF采用运动波(Kinematic Wave)近似表达均质非饱和带中的一维理查德方程,这种描述方程更适用与大空间尺度模型(Niswonger 等,2006)。

Brooks-Corey 模型(Brooks和Corey,1964)用来描述非饱和带渗透系数和土壤含水量之间的关系:

式中Ks为饱和渗透系数;θr为残余含水量;θs为饱和含水量;ε为Brooks-Corey 常数。   本次研究的主要目的是估计区域地下水补给,此过程中将非饱和带在垂向上假设成均质条件一般认为是合理的(Niswonger和Prudic,2009)。UZF的主要输入参数包括入渗强度(FINF),垂向渗透系数(VKS),饱和含水量(THTS) 和Brooks-Corey系数(EPS)(见表5-1). 上述土壤水均衡计算得到的正值Rp作为入渗量输入,而负值Rp表示水分从根植层下部非饱和带向根植层运动,作为蒸散发输入模型。 

5.3.3 模型数据输入

蒸发量采用第5章中CRAE模型计算的实际蒸发量。土壤质地(包括砂,壤土和粘土含量)及土壤容重从中国1比100万土壤数据库获得(图 5-1)。该数据库为分辨率1km的栅格结构,按照模型网格内主要质地类型提升到模型适用2km网格。van Genuchten模型(van Genuchten,1980)中的七个土壤水力参数利用Schaap(2001)利用美国土壤数据开发的土壤转换函数计算得到。θs和Ks(图 5-2)对应与UZF中的VKS和THTS。Brooks-Corey 常数根据Morel-Seytoux 等(1996)导出公式计算:

土壤含水量数据从中国农业气象数据库获得,该数据库在华北平原范围内包括104个观测站,观测数据为10,20和50 cm处的土壤重量含水量。每月8号、18号和28号进行观测。数据为相对含水量(土壤含水量与田间持水量之比),利用土壤容重和田间持水量数据将相对重量含水量转换为体积含水量后输入模型。  华北平原没有详细精确灌溉量分布可供利用,小麦-玉米轮作灌溉需水量约为4,200 m3/ha/yr(柳长顺等,2004)。华北平原小麦一般灌溉4次或5次,玉米一般灌溉3或4次。模型将灌溉处理成7次,每年3月至7月各一次,9月和10月各一次,不区分作物类型,每次灌溉量为60 mm。 

5.3.4 与饱和带水流模型耦合

根据上述讨论,本次研究中非饱和带-饱和带耦合水流模型实际上可认为是一个三层模型。第一层模型为根植层的土壤水均衡模型,该模型将土壤水储变量和蒸散发从降水量中扣除,剩余量作为根植层底界和下部非饱和带的水分交换量;第二层模型为非饱和水流数值模型,根植层底界渗出水分通过该模型传输到下部饱和带地下水;第三层为饱和带水流模型,其表达与第5章中区域地下水模型一致。 模型模拟时段为1993年1月至2008年12月,每月作为1个应力期。模型网格剖分与区域地下水模型一致。模型初始条件采用区域非稳定水流模型计算得到的1992年底的水位分布。非稳定流模型中开采量为年地下水总开采量,因耦合模型认为农业灌溉开采量占总开采的70%,并平均分配到3月至7月,9月及10月。其余30%开采量作为工业开采和生活开采量平均分配到各月。地下水开采的空间分布与非稳定流模型一致,而灌溉量均匀分布输入模型。 非饱和带中水流过程模拟中,初始含水量的确定通常采用非饱和带稳定流模拟,或者反复运行模型较长时间获得(Keese等,2005;Wang等,2009)。然而,由于华北地下水已经处于超采状态,不可能建立现状开采条件下的稳定流模型。本次研究中,采用0.2为模型初始含水量,从1960开始,运行模型至1992年,得到的含水量分布作为研究模拟期的初始含水量。因为开采初期华北平原地下水埋深范围1-2 m,1960年-1992年模拟期间由于初始含水量带来的计算误差将会很快被修正。 

5.4 模拟结果分析

5.4.1 计算区域地下水补给量

计算年平均地下水补给量为153 mm,占年降水和灌溉水量的18%。华北平原范围内计算地下水补给量在0-260 mm/yr 范围内的面积占整个平原面积的90%。地下水补给的时空分布与降水量、灌溉量和蒸散量的时空分布高度相关。较强地下水补给发生在7月至9月(图 5-3)。图 5-3表明初始含水量,降水、灌溉和蒸散发对地下水补给有显著的季节效应。冬季地下水补给量与春季相比,春季虽然降水量较大而且是主要灌溉时期,但春季蒸发量也较大,最终产生的地下水补给量比冬季反而要小。年均地下水补给量在山前平原地区较高,为0-360 mm/yr(图 5-4),占年降水量和灌溉量的44%。地下水补给量低于平原均值的地区约占平原面积的60%,大部分位于中部平原和滨海平原区。因此,山前平原相对较高的地下水补给使整个平原地下水补给的平均值较中东部平原的真实值偏高。计算结果显示的地下水补给量的空间分布模式与前人研究成果(如,von Rohden 等,2010;Foster 等,2004)是一致的。年均净补给量(补给量与地下水蒸发量之差,地下水蒸发截止深度按4 m计算)为117 mm。 

土壤质地类型是影响地下水补给量的主要因素,从图 5-5可以看出,粗粒径土壤地区的地下水补给量明显高于其他地区:壤砂土分布地区的地下水补给量为~9 - ~360 mm/a(分别占占年降水和灌溉量之和的~1%和~ 43%);壤土分布地区地下水补给量为~1 - ~270 mm/yr(分别占占年降水和灌溉量之和的~ 0.1% 和~ 32 %);粘性土分布地区的地下水补给量为~ 0.1 - ~220 mm/yr,(分别占年降水和灌溉量之和的约0.01%和约26%)。由于降水量和灌溉量的空间分布不均,计算地下水补给量壤土和粘土,壤砂土和壤土分布区之间都存在重叠。同一土壤质地分布地区内的地下水补给的分布离散范围都大于不同质地之间的离散范围,说明虽然土壤质地是影响地下水补给量的重要要因素,却不是造成地下水补给空间变异性的控制因素。 

地下水入渗补给系数(Rc)定义为计算地下水入渗补给与当年降水量与灌溉量加和之比,Rc=R/(P+I),用来定量描述补给量与非饱和带厚度变化的相关关系。图 5-6显示了各土壤质地类型分布地区不同非饱和带增厚范围内的入渗补给系数。可以看出,非饱和带增厚对地下水补给的影响应与非饱和带的渗透性有关。对于壤砂土地区,非饱和带厚度小于30 m时,渗透系数没有呈现随非饱和带增厚而逐渐减小的现象(图 5-6a)。主要原因是对于深部非饱和带,含水量主要与质地类型有关,壤砂土一般含水量较小,非饱和带储存量变化对地下水补给(最终达到潜水位的入渗量)影响不大。对于壤土和粘土,渗透系数随非饱和带厚度的增加而减小(图 5-6b,c)。渗透系数和非饱和带厚度之间的相关关系在前人的场地试验研究中已有研究,比较一致认为随非饱和带厚度增加,渗透系数存在一先增大后减小至稳定的过程(如,李亚峰等,2007;李金柱,2009)。对于整个华北平原,非饱和带厚度对模拟内地下水入渗补给系数大小影响不大(图 5-6)。因为入渗补给系数的定义是补给量和当年的降雨量和灌溉量的比值,不能够反映非饱和带增厚对地下水补给的缓冲和滞后作用,因而不能直观显示非饱和带增厚对地下水补给的影响。

从模型中选取分别位于山前平原和中部平原的四个模型单元(见图 5-4),涵盖了华北平原主要土壤质地类型,山前平原的三个模型单位均为大厚度(> 30 m)非饱和带,分别统计其地下水补给量随时间变化(图 5-7)。模拟期间,山前平原三个观测点非饱和厚度从约20 m增厚至近50 m,中部平原观测点非饱和带厚度变化范围为1.5 - 3 m。从图可以看出,非饱和厚度对地下水补给量在时间域内的分布和大小都有显著影响。大厚度非饱和带(地下水位埋深大)不但使地下水补给到达时间发生滞后,并对其时间域内分布造成平滑现象。由于流经非饱和带时部分水分转换成非饱和水分储存量(土壤含水量增大),最终到达潜水面的地下水补给量也相应有所减少。然而,降水量和灌溉变化造成的非饱和带储存量变化要大于包气带增厚造成的储存量变化,例如,对于三个大厚度包气带模型单元1996年降水量增大造成的非饱和储存量增加量占到了整个模拟期间非饱和带储存量增加量的30-80%。 图 5-7也显示,对于大厚度包气带,当包气带上部入渗水量没有剧烈变化时,地下水补给量可以认为是一定值:2004年至2008年之间,沙壤土模型单元年地下水补给量稳定在260-300 mm之间,均值约为280 mm(图 5-7a)。此时,非饱和带内土壤水储变量可以忽略不计,地下水补给量可以用根植层土壤水均衡计算得到的潜在地下水补给表示。该结果与前人的试验研究一致(Wu等,1996)。总之,对于大厚度包气带地区,当非饱和带的质地结构和含水量饱和情况不明时,直接利用土壤水均衡计算的地下水入渗量作为地下水补给量将会与实际情况产生较大偏差,非饱和带储存量变化也是不容忽视的重要因素。 

5.4.2 模型参数灵敏度分析

本小节分析地下水补给系数对模型各参数的灵敏度,确定对地下水补给量影响最大的模型参数。模型灵敏度分析通过对模型初始含水量和各模型参数在一定范围内变化(其他参数不变)的情况下模型结果的变化情况,反映最为敏感的参数(图 5-8)。具体为对饱和含水量(THTS)分别乘以(0.1,0.2,0.5,1.5,2.0,2.5)内各乘数,对初始土壤含水量(THTI)分别乘以(0.5,0.7,0.8,1.5,2.0,2.5)内各乘数;对给水度(Sy)和垂向饱和渗透系数(VKS)分别乘以(0.1,0.2,0.5,2.0,5.0,10.0)内各乘数,对Brooks-Corey常数分别乘以(0.3,0.4,0.5,2.0,5.0,10.0)内不同乘数,某一参数变化时其他参数保持不变。Brooks-Corey 常数乘以0.3和0.4时,模型未能收敛,因此没有参与分析。结果显示,计算地下水补给量对初始土壤含水量最为敏感,将初始含水量增大2.5倍时将使地下水补给量增大10倍。说明土壤初始含水量对地下水补给的重要性,初始含水量较低时,入渗水分相当一部分将转换成非饱和带水分储存量,因而降低了地下水补给,反之,较高的初始含水量将造成较大的地下水补给。这也说明,因为降水量和灌溉的时间分布必然影响非饱和带的前期含水量,因此也将决定地下水补给过程和补给量大小。对补给量非常敏感的另外一个参数是给水度。这是因为UZF不仅利用给水度计算饱和带地下水储存量的变化,也可利用其计算湿润锋面的运动速度和潜水面上升时的储存量增加。地下水补给与Brooks-Corey常数成负相关关系,当参数值减小时,Brooks-Corey常数比饱和含水量和饱和渗透系数对地下水补给量更为敏感。  

5.4.3 非饱和带地下水滞留时间

包气带中地下水运动时间反映了地下水对地表条件(如气候条件改变,补给条件改变和污染物释放)改变的响应快慢程度(Flint等,2004;Grismer,2013)。包气带中地下水运动时间受控于地下水渗漏量,包气带厚度和岩性(特别是渗流途径上低渗透性介质)和有效孔隙度(Scanlon等,2006),包气带中地下水运动时间可按照下式进行计算(Nichols和Freshley,1993):

式中t为包气带中地下水运动时间;ne为有效孔隙度;Θ为有效饱和度;Z为包气带厚度;q为垂向达西流速,也即地下水补给率。有效饱和度(即无量纲化的含水量)定义为:

可通过上述非饱和带水流模型计算得到q和Θ,而有效孔隙度通常难以获得,计算过程中根据包气带岩性采用经验值。本文计算过程中有效孔隙度采用0.25。

计算非饱和带中地下水运动时间的空间分布见图 5-9。区域尺度上,利用上述方法计算得到的非饱和带地下水运动时间在量级上与碳-14和氚同位素反应的时间尺度是一致的。华北平原浅层地下水中14C含量(现代碳百分比)一般小于50pmc(刘存富等,1997),深度小于50 m的浅井中地下水中氚含量一般大于15TU(Chen等,2005)。 

图 5-10显示了沿石家庄-衡水剖面计算包气带地下水运动时间和非饱和带厚度。然而,山前地区非饱和带地下水运动时间与前人研究的3H-3He地下水年龄(von Rohden等,2010)存在较大不一致。石家庄地区地下水漏斗附近地下水年龄明显小于包气带中的地下水运动时间。然而,该地区地下水年龄随水流路径变老的现象也不明显,合理的解释应与山前地区存在的山区地下水侧向补给有关,地下水开采在加速侧向补给的同时,也使垂向补给的地下水与侧向补给的地下水的混合加剧。靠近中部平原地区的浅层地下水年龄更接近于包气带中的地下水运动时间,说明该地区地下水补给的主要来源为垂向补给。

5.4.4 模拟结果验证

图 5-11显示了105眼观测井处共13900个观测点的观测水位与计算水位对比。两者间均方根误差 (RMSE) 约为10 m,考虑到模拟时段内的地下水动态变化范围(-18.2 - 92.5 m),拟合精度是可以接受的。模拟时段内华北平原平均降水量约为520 mm,对于此类干旱地区,117 mm/yr的净地下水补给量对于通常认识(如 Scanlon 等,2006)来说是过高的。前人一般认为华北平原降水入渗形成的地下水补给占降水量的~1 % - ~20 % (陈望和,1999),这说明部分地下水补给量应该来自于灌溉回归补给(Kendy等,2004)。地下水开采量(~160mm/yr)与地下水净补给之差约为40 mm/yr,面积加权平均给水度按0.075计算,平均潜水位下降速度为0.5 m/yr,与水位观测结果一致。

5.5 本章小结

采用土壤水均衡模型和非饱和带水流数值模拟的方法计算了华北平原地下水入渗补给量的时空变化特征。该方法的基本原理是利用根植层的土壤水均衡计算根植层底部水分通量,并将其通过一维非饱和带水流模型传送至潜水面,饱和带水流模型为非饱和带水流运动模拟提供下边界条件。既避免了复杂的大气-植被-土壤水复杂的模拟过程,也能够快速有效的直接模拟大厚度非饱和带地区非饱和带对地下水补给过程的影响。华北平均计算地下水补给量在山前平原地区最大约为360 mm/yr,向中东部平原降低,最高为至 260 mm/yr。整个平原区年均地下水补给量为153 mm。计算结果与前人利用不同方法(如环境示踪剂方法)计算的局部地下水补给范围一致。模拟结果表明,华北平原非饱和带增厚既对地下水入渗补给起到显著的缓冲和滞后作用,由于厚度增加带来的非饱和带储存量增大也一定程度上减小了地下水补给。而利用入渗补给系数的概念则不能反映非这些饱和带对地下水补给的影响作用。 模型结果显示大厚度非饱和带对地下水补给时间域内的分布及大小的影响因土壤质地而异:非饱和带厚度的增大对壤砂土地区的地下水补给影响不显著,但会明显减小壤土和粘土地区的地下水补给量。模拟结果表明,土壤质地类型和非饱和带厚度对华北平原地下水补给的控制作用要低于降水量和灌溉的时空分布。 本研究中采用的土壤水均衡方法,而没有利用非饱和带水流模型来计算蒸发量,因此也避免了对计算蒸发量进行校正,也因此避免了传统非饱和带水流模型需要的植被信息,也为利用其他更适用于区域蒸发量的计算方法(如遥感方法)提供了模型接口,对于其他类似地区进行快速高效地计算区域地下水补给分布,建立耦合非饱和带-饱和带水流模型有重要的借鉴意义。