北京城近郊区地下水数值模拟模型及应用

北京城近郊区地下水数值模拟模型及应用北京城近郊区地下水数值模拟模型及应用王新娟1,2崔亚莉2邵景力21.北京市地质工程勘察院水资源研究所1000372.中国地质大学（北京） 100083摘要：本文在深入了解北京城近郊区地下含水层及其岩性的基础上建立水文地质概念模型，运用FEFLOW软件建立了研究区三维地下水流数值模拟模型，并对模型进行了识别和检验。通过模型进行地下水补给资源和可开采资源评价，提出地下水可开采资源量。关键词：模拟模型水源地初始流场模拟流场引言北京市作为我国的政治、经济、文化中心，同时又是一个严重缺水的城市。由于长期以来缺1王新娟（1973-），女，工程师，中国地质大学（北京）水资源与水文学专业博士,现从事水工环评价与模型研究。乏科学的管理、盲目过量集中开采地下水，导致城近郊区地下水位出现区域性大幅度下降，部分地区含水层疏干或半疏干，导致了地下水质恶化，水源污染，地面沉降，地下水自然升温等各种水文、环境、工程地质问题，加剧了水资源供需矛盾，目前水资源已成为制约首都经济发展、城市功能正常发挥的重要因素之一。因此，在北京地区建立地下水资源可持续开采模型，通过调参、计算解决上述问题势在必行。1水文地质条件永定河冲洪积扇地区地处华北平原的西北端，地势由西北向东南倾斜，地面高程100〜45米，冲洪积扇顶部地形坡度3%。，至南部坡度为1%。。第四系沉积规律主要受永定河河道变化的影响，八宝山以北地区，第四系沉积厚度可达200米以上。含水层的主要作用是接受大气降水、河谷潜流水、山前侧向迳流水、河渠、灌溉及人工回灌水的入渗补给，水在含水层中传输运移，经历许多水文循环过程。植被吸收、蒸发、包气带消耗、人工开采等在系统中均可产生不同作用。含水层的岩性和结构可将系统划分成三个子系统（1）单层砂卵石储水区，位于莲花池-昆明湖以西；含水层渗透性能好，渗透系数300—500m/d。（2）二到三层含卵石砂砾石储水区，主要位于门头村一八里庄一陶然亭一马家堡一带，渗透系数50—300m/do（3）多层中细砂及含砾砂储水区，该区分布在上述地区以东的广大地区，渗透系数30—50m/do2水文地质概念模型和数值模型2.1含水层结构及概化为了建立城近郊区地下水流模型，首先要对实际水文地质条件概化，建立水文地质概念模型。根据前述水文地质条件，模型在空间上分为三层，第一含水层组（主要是部分农业开采，地面以下20m左右）、弱透水层（地面以下20〜40m左右的一个粘性土层）、第二含水层组（主要是城市工业和生活用水及部分农业开采，地面以下40〜130m左右）。研究区含水层系统涉及的边界主要是流量边界。就整体而言，西北部边界是侧向补给边界，主要接受山前侧向补给；地下水初始流场如图1—1示：研究区的东南部边界为整个含水层系统的天然排泄边界，垂向边界包括永定河河床潜流及河道渗漏补给，永定河引水渠和京密引水渠大部分有衬砌，仅有少量渠道水入渗地下；大气降水入渗补给是系统的主要补给资源之一，地下水的蒸发主要指潜水的自然蒸发，据北京市水文公司试验资料：潜水埋深大于3.5m时，各类岩性的潜水蒸发折算系数趋近于零；底部边界和顶托补排，系统的第四系含水层底部，大多数为透水性差的第三系地层，预计其对系统内部影响很小，但在海淀部分地区，奥陶系石炭系二迭系侏罗系地层与第四季直接接触，地下水垂向顶托补排作用肯定存在；农田灌溉用水包括地下水和地表水，灌溉水除蒸发、植物吸收、包气带损耗外，还有一部分水入渗补给地下水。综上所述，研究区可概化成非均质各向同性、空间三维结构、非稳定地下水流系统，即地下水系统的概念模型。2.2地下水数值模型对于上述非均质、各向同性、空间三维结构、非稳定地下水流系统，可用如下微分方程的定解问题来描述：式中：Q—分别为含水层渗流区域；Q—渗流区域；h—h=h（x,y,z），含水层的水位标高（m）;Kx、Ky、Kz—分别为x、y、z方向的渗透系数（m/d）;Kn—边界面法向方向的渗透系数（m/d）;S—自由面以下含水层储水系数（1/m）;卩一潜水含水层在潜水面上的重力给水度；（T—河流底部弱透水层的阻力系数，（T=L/KS,L为底部弱透水层的厚度，Ks为河流底部弱透水层垂向渗透系数（m/d）;—含水层的源汇项（1/d）;p—潜水面的蒸发和降水等（1/d）;ho—含水层的初始水位分布（m），ho=ho（x,y,z）;r—渗流区域的上边界，即地下水的自由表面；r—渗流区域的下边界，即承压含水层底部的隔水边界；r—渗流区域的侧向边界；暄一边界面的法线方向；q（x,y,z,t）—定义为二类边界的单宽流量（m2/d・m），流入为正，流出为负，隔水边界为 0;研究区进行网格剖分，剖分后的模拟区共17000个结点，4327个单元格。2.3模型的检验和识别模型的识别与检验过程是整个模拟中极为重要的一步工作，通常要进行反复地修改参数和调整某些源汇项才能达到较为理想的拟合结果。此模型的识别与检验过程采用的方法也称试估一校正法，它属于反求参数的间接方法之一。通过拟合同时期的流场和长观孔的历时曲线， 识别水文地质参数、边界值和其它均衡项，使建立的模型更加符合研究区的水文地质条件，以便更精确地定量研究模拟区的补给与排泄，预报给定水资源开发利用方案下的地下水位。对研究区地下水系统进行了识别和验证。通过反复调整参数和均衡量，识别水文地质条件，确定了模型结构、参数和均衡要素。模拟末期（1995年12月31日）含水层的模拟流畅与实际流畅对比见图1—2;长观孔的历时曲线模拟见图1—3（略）通过1995年模拟模型的研究，建立模型，对模型进行验证和检验，模型拟合效果较好，为了进一步验证模拟模型的可行性，利用2000年1月1日-12月31日一年的实测资料，分12个时段进行验证，验证结果：2000年是枯水年，地下水位下降很大，验证期末的潜水地下水模拟流场与实际流场对比见图1-4。（略）由图示可见，所建立的模拟模型基本达到模型精度要求，符合水文地质条件，基本反映了地下水系统的水力特征，可利用模型进行地下水位预报2.4地下水模拟模型的应用应用已建立的地下水流模型，可进行研究区地下水资源量评价、提出北京地区地下水可持续利用条件下的开采量。研究区模拟期内地下水的补给来源述要有三个方面，即大气降水、山前侧向径流、地表水转化量（包括永定河弃水、引水渠渗漏），其中降水入渗量占31.90%，各种地表水转化量占29%左右，其中永定河弃水入渗量占21.53%。1995年总补给量为64384.2310忖，按研究区面积1082.89km2计算，其平均补给模数为59.5104m3/km2?a。1995年降水量为平偏丰水年，加上永定河弃水补给地下水，整个地下水系统呈现正均衡，盈1257.41104m3。为此，有必要利用多年降水量资料，开展不同保证率下的地下水资源量评价。2.4.1潜水地下水资源量评价（1）补给资源量由地下水均衡分析结果可知，潜层地下水补给来源于降水入渗、河渠入渗、灌溉回渗、以及侧向边界补给。其中占主导地位的是降水入渗、河渠入渗、侧向补给以及灌溉回渗，合计约占总补给量的87.4%。选取石景山气象站资料降水量观测数据，按保证率计算公式p，给出1977-2000年降水量保证率（表4—1），公式中：P—经验保证率（%）;m—年降水量按大小排列序列号；n—统计年数。由表可见，能代表50%降雨保证率的典型年为1991年、1995年，降水量平均为625.92mm,由此可以推知其余地区在同一保证率下的降水量，这样可以计算出研究区在该保证率下总的降雨入渗量，为20538.87104m3。侧向补给等各项补给量见表4-2,在有永定河弃水的条件下，系统总补给量为60132.23104m3，平均补给模数为55.58104m3/km2a,若不考虑永定河弃水，则总补给量为46271.97104m3。能代表75%保证率的典型年为2000年，降水量为428.6mm，由此可推知相应地区在同一保证率下的降水量，从而计算出研究区在 75%保证率下总的降水入渗量为10443.2104m3。侧向补给量为6435.7104m3,河渠入渗量、灌溉回渗根据2000年模型验证，计算结果分别为3969.2104m3、2931.80104m3,见表4-2。（略）表4-1 石景山气象站降水量保证率计算表序号年份降雨量（mm）频率（％）序号年份降雨量（mm）频率（%）11977877.14.00131995539.052.0021994775.28.00141989486.256.0031985752.412.00151992482.060.0041986728.216.00161983476.964.0051996689.720.00171982461.668.00

补给项平水年(50%保证率)枯水年(75%保证率)补给量(10W)比例补给量(104m3)比例降水入渗量20538.8734.1610443.238.90永定河弃水入渗量13860.2623.056435.723.97引水渠入渗量4748.007.902931.8010.92灌溉回渗量6871.1111.432405.008.96山区流入量10248.7817.041564.25.83灰岩顶托补给量3865.216.433865.2111.42合计60132.2310027645.11100.00面积10821082补给模数55.5825.5561998689.724.00181993447.372.007198761998689.724.00181993447.372.0071987666.428.00192000428.676.0081988659.532.00201981392.080.0091978655.436.00211999375.484.00101979631.240.00221984358.188.00111990613.844.00231997351.892.00121991609.848.00241980271.596.00平均9.12 55丿、i 可知―平Si证均m/aMk量计算表-21水东水朝阳。潜水地下水资源补给(2)储存量根据潜水含水层给水度分区，分别计算出各区的枯水位以下含水层的平均厚度和分布面积，采用下式计算潜水储存量：Q储=F?H?□式中：Q储一含水层的储存量（104m3）;F—计算区面积（104m2）;H—枯水位以下含水砂层平均厚度（m）;□—潜水含水层的给水度。计算得到研究区内潜水含水层储存量为57.461108m3（表4-3）,由此可见，潜水含水层组具有较强的可调节性。表4—3 研究区潜水各个分区的储存量计算表ID面积（km2）含水层厚度（m）给水度存储量(108m3)13.3813.560.010.005292.9350.770.2310.852387.5350.770.2310.221419.3839.430.261.987543.6350.970.163.558613.5025.970.1860.6527172.0120.110.237.9568378.6019.930.1410.564934.2520.830.251.7841023.5418.300.220.9481156.9419.880.2152.43412163.8220.890.196.50257.4612.4.2承压地下水资源评价承压水位于研究区中东部的含水层，受气候的影响相对潜水含水层小，所以补排量主要是人工开采、层间排泄、山区单一含水层的越流补给和迳流。研究区内迳流非常较为缓慢，迳流补给和排泄，人工开采和侧向流出是地下水主要的排泄方式，而越流补给和层间补给是承压水的主要补给来源。1995年承压水含水层组总补给量为5414.20104m3,总排泄量为59602.2104m3,均衡差为-54188.00104m3。深层水补给量主要是浅层越流和迳流补给，弹性释水量计算公式为：Q弹=厶h?F?lxc式中：Q弹一弹性释水量（104m3）;△h—静水位与含水层顶板之差（m）;F—计算区面积（104m2）;xc—含水层弹性释水系数（无量纲）经计算得到研究区内深层水的弹性释水量为89660.00104m3（表4-4）。

表4-4 研究区深层各个分区储存量计算表编号面积（km2）含水层厚度（m）释水系数存储量（108m3）13.380.000.000.00292.930.000.000.00387.530.000.000.00456.9440.780.02010.4675163.8283.420.0192.597619.380.000.000.00743.620.000.000.00813.583.190.0220.247923.5463.490.01860.27711378.6045.400.0142.40612172.0165.690.0232.5991334.2543.470.0250.372合 计8.9662.4.3可采资源量地下水可开采资源量是指在经济和技术合理的条件下，开采过程中不发生水质恶化或其它不良地质现象（如地面沉降、水的硬度升高等），并对生态环境不致造成不利影响的情况下，有保证的地下水资源量。研究区西部潜水含水层组较厚，约20〜40m，含水层渗透性很强，是一个巨大的地下水库，具有良好的多年调节作用；因此，可以将平水年地下水补给量可开采资源量计算的基础。1995年研究区开采量为6.06108m3，平均补给量为6.01108m3,全区基本平衡。但是应当注意：西郊地区即是地下水的补给区，也是北京市水源地和工业自备井集中开采区，长期的过量开采，已导致局部地区地下水位大幅度下降，含水层疏干状态。虽然西郊局部地区存在地下水超采问题，但是大部分地区具有很好的调蓄功能， 可以利用永定河、砂石坑等进行人工回灌，因此确定地下水可开采资源量为5.5108m3/a是有保证的。3结语本研究在水文地质条件分析的基础上，运用FEFLOW建立了北京城近郊区地下水数值模拟模型，通过模型模拟调参，对模型进行识别和验证。计算出研究区平水年地下水补给资源和可开采资源评价分别为6.013X108m3/a和5.5X108m3/a。北京城近郊区地下水数值模拟模型及应用北京城近郊区地下水数值模拟模型及应用参考文献：吕晓检，《北京市城近郊区地下水管理模型》1991.8Bear,J.,许涓铭等译，地下水水力学，北京,1985.许恒力等，水资源开发与保护，地质出版社,北京，2001.宾德智等，北京市中心区供水水文地质勘察研究报告，1964.8Anderson,M.P.andW.W.Woessner,AppliedGroundwaterModelingSimulationofFlowandAdvectiveTransport,AcademicPress,1992.王新娟1,2崔亚莉2邵景力21.北京市地质工程勘察院水资源研究所 1000372.中国地质大学（北京） 100083摘要：本文在深入了解北京城近郊区地下含水层及其岩性的基础上建立水文地质概念模型，运用FEFLOW软件建立了研究区三维地下水流数值模拟模型，并对模型进行了识别和检验。通过模型进行地下水补给资源和可开采资源评价，提出地下水可开采资源量。关键词：模拟模型水源地初始流场模拟流场引言北京市作为我国的政治、经济、文化中心，同时又是一个严重缺水的城市。由于长期以来缺乏科学的管理、盲目过量集中开采地下水，导致城近郊区地下水位出现区域性大幅度下降，部分地区含水层疏干或半疏干，导致了地下水质恶化，水源污染，地面沉降，地下水自然升温等各种水文、环境、工程地质问题，加剧了水资源供需矛盾，目前水资源已成为制约首都经济发展、城市功能正常发挥的重要因素之一。因此，在北京地区建立地下水资源可持续开采模型，通过调参、计算解决上述问题势在必行。1水文地质条件永定河冲洪积扇地区地处华北平原的西北端，地势由西北向东南倾斜，地面高程100〜45米，冲洪积扇顶部地形坡度3%。，至南部坡度为1%。。第四系沉积规律主要受永定河河道变化的影响，八宝山以北地区，第四系沉积厚度可达200米以上。含水层的主要作用是接受大气降水、河谷潜流水、山前侧向迳流水、河渠、灌溉及人工回灌水的入渗补给，水在含水层中传输运移，经历许多水文循环过程。植被吸收、蒸发、包气带消耗、人工开采等在系统中均可产生不同作用。含水层的岩性和结构可将系统划分成三个子系统（1）单层砂卵石储水区，位于莲花池一昆明湖以西；含水层渗透性能好，渗透系数300—500m/d。（2）二到三层含卵石砂砾石储水区，主要位于门头村一八里庄一陶然亭一马家堡一带，渗透系数50—300m/do（3）多层中细砂及含砾砂储水区，该区分布在上述地区以东的广大地区，渗透系数30—50m/d。2水文地质概念模型和数值模型2.1含水层结构及概化为了建立城近郊区地下水流模型，首先要对实际水文地质条件概化，建立水文地质概念模型。根据前述水文地质条件，模型在空间上分为三层，第一含水层组（主要是部分农业开采，地面以下20m左右）、弱透水层（地面以下20〜40m左右的一个粘性土层）、第二含水层组（主要是城市工业和生活用水及部分农业开采，地面以下40〜130m左右）。研究区含水层系统涉及的边界主要是流量边界。就整体而言，西北部边界是侧向补给边界，主要接受山前侧向补给；地下水初始流场如图1—1示：研究区的东南部边界为整个含水层系统的天然排泄边界，垂向边界包括永定河河床潜流及河道渗漏补给，永定河引水渠和京密引水渠大部分有衬砌，仅有少量渠道水入渗地下；大气降水入渗补给是系统的主要补给资源之一，地下水的蒸发主要指潜水的自然蒸发，据北京市水文公司试验资料：潜水埋深大于3.5m时，各类岩性的潜水蒸发折算系数趋近于零；底部边界和顶托补排，系统的第四系含水层底部，大多数为透水性差的第三系地层，预计其对系统内部影响很小，但在海淀部分地区，奥陶系石炭系二迭系侏罗系地层与第四季直接接触，地下水垂向顶托补排作用肯定存在；农田灌溉用水包括地下水和地表水，灌溉水除蒸发、植物吸收、包气带损耗外，还有一部分水入渗补给地下水。综上所述，研究区可概化成非均质各向同性、空间三维结构、非稳定地下水流系统，即地下水系统的概念模型。2.2地下水数值模型对于上述非均质、各向同性、空间三维结构、非稳定地下水流系统，可用如下微分方程的定解问题来描述：0式中：Q—分别为含水层渗流区域；Q—渗流区域；h—h=h（x,y,z），含水层的水位标高（m）;Kx、Ky、Kz—分别为x、y、z方向的渗透系数（m/d）;Kn—边界面法向方向的渗透系数（m/d）;S—自由面以下含水层储水系数（1/m）;卩一潜水含水层在潜水面上的重力给水度；（T—河流底部弱透水层的阻力系数，（T=L/Ks,L为底部弱透水层的厚度，Ks为河流底部弱透水层垂向渗透系数（m/d）;—含水层的源汇项（1/d）;p—潜水面的蒸发和降水等（1/d）;ho—含水层的初始水位分布（m）,ho=ho（x,y,z）;r—渗流区域的上边界，即地下水的自由表面；r—渗流区域的下边界，即承压含水层底部的隔水边界；r—渗流区域的侧向边界； 萨边界面的法线方向;q（x,y,乙t）—定义为二类边界的单宽流量（m2/d.m），流入为正，流出为负，隔水边界为 0;研究区进行网格剖分，剖分后的模拟区共 17000个结点，4327个单元格。2.3模型的检验和识别模型的识别与检验过程是整个模拟中极为重要的一步工作，通常要进行反复地修改参数和调整某些源汇项才能达到较为理想的拟合结果。此模型的识别与检验过程采用的方法也称试估一校正法，它属于反求参数的间接方法之一。通过拟合同时期的流场和长观孔的历时曲线，识别水文地质参数、边界值和其它均衡项，使建立的模型更加符合研究区的水文地质条件，以便更精确地定量研究模拟区的补给与排泄，预报给定水资源开发利用方案下的地下水位。对研究区地下水系统进行了识别和验证。通过反复调整参数和均衡量，识别水文地质条件，确定了模型结构、参数和均衡要素。模拟末期（1995年12月31日）含水层的模拟流畅与实际流畅对比见图1—2;长观孔的历时曲线模拟见图1—3（略）通过1995年模拟模型的研究，建立模型，对模型进行验证和检验，模型拟合效果较好，为了进一步验证模拟模型的可行性，利用2000年1月1日—12月31日一年的实测资料，分12个时段进行验证，验证结果：2000年是枯水年，地下水位下降很大，验证期末的潜水地下水模拟流场与实际流场对比见图1-4。（略）由图示可见，所建立的模拟模型基本达到模型精度要求，符合水文地质条件，基本反映了地下水系统的水力特征，可利用模型进行地下水位预报2.4地下水模拟模型的应用应用已建立的地下水流模型，可进行研究区地下水资源量评价、提出北京地区地下水可持续利用条件下的开采量。研究区模拟期内地下水的补给来源述要有三个方面，即大气降水、山前侧向径流、地表水转化量（包括永定河弃水、引水渠渗漏），其中降水入渗量占31.90%，各种地表水转化量占29%左右，其中永定河弃水入渗量占21.53%。1995年总补给量为64384.2310务3，按研究区面积1082.89km2计算，其平均补给模数为59.5104m3/km2?a。1995年降水量为平偏丰水年，加上永定河弃水补给地下水， 整个地下水系统呈现正均衡，盈1257.41104m3。为此，有必要利用多年降水量资料，开展不同保证率下的地下水资源量评价。2.4.1潜水地下水资源量评价（1）补给资源量由地下水均衡分析结果可知，潜层地下水补给来源于降水入渗、河渠入渗、灌溉回渗、以及侧向边界补给。其中占主导地位的是降水入渗、河渠入渗、侧向补给以及灌溉回渗，合计约占总补给量的87.4%。选取石景山气象站资料降水量观测数据，按保证率计算公式兰二［，给出1977-2000年降水量保证率（表4—1），公式中：P—经验保证率（%）;m—年降水量按大小排列序列号；n—统计年数。由表可见，能代表50%降雨保证率的典型年为1991年、1995年，降水量平均为625.92mm,由此可以推知其余地区在同一保证率下的降水量，这样可以计算出研究区在该保证率下总的降雨入渗量，为20538.87104m3。侧向补给等各项补给量见表4-2,在有永定河弃水的条件下，系统总补给量为60132.23104m3，平均补给模数为55.58104m3/km2a,若不考虑永定河弃水，则总补给量为46271.97104m3。能代表75%保证率的典型年为2000年，降水量为428.6mm，由此可推知相应地区在同一保证率下的降水量，从而计算出研究区在 75%保证率下总的降水入渗量为10443.2104m3。侧向补给量为6435.7104m3,河渠入渗量、灌溉回渗根据2000年模型验证，计算结果分别为3969.2104m3、2931.80104m3,见表4—2。(略)表4-1 石景山气象站降水量保证率计算表序号年份降雨量(mm)频率(%)序号年份降雨量(mm)频率(%)11977877.14.00131995539.052.0021994775.28.00141989486.256.00

31985752.412.00151992482.060.0041986728.216.00161983476.964.0051996689.720.00171982461.668.0061998689.724.00181993447.372.0071987666.428.00192000428.676.0081988659.532.00201981392.080.0091978655.436.00211999375.484.00101979631.240.00221984358.188.00111990613.844.00231997351.892.00121991609.848.00241980271.596.00均 平9.12 55FT在P，年-,COL,T多0水。情FT在P，年-,COL,T多0水。情牽1分0上证「令O以魁?a00陳a24下O2/IKV在？i率？,量55^m知SI明-4證4密!-究均弊补给项平水年（50%保证率）枯水年（75%保证率）补给量（10m3）比例补给量(104m3)比例降水入渗量20538.8734.1610443.238.90永定河弃水入渗量13860.2623.056435.723.97引水渠入渗量4748.007.902931.8010.92灌溉回渗量6871.1111.432405.008.96山区流入量10248.7817.041564.25.83灰岩顶托补给量3865.216.433865.2111.42合计60132.2310027645.11100.00面积10821082补给模数55.5825.55（2）储存量根据潜水含水层给水度分区，分别计算出各区的枯水位以下含水层的平均厚度和分布面积，采用下式计算潜水储存量：Q储=F?H?□式中：Q储一含水层的储存量（104m3）;F—计算区面积（104m2）;H—枯水位以下含水砂层平均厚度（m）;□—潜水含水层的给水度。计算得到研究区内潜水含水层储存量为57.461108m3（表4-3）,由此可见，潜水含水层组具有较强的可调节性。表4—3 研究区潜水各个分区的储存量计算表ID面积（km2）含水层厚度（m）给水度存储量（108m3）13.3813.560.010.005292.9350.770.2310.852387.5350.770.2310.221419.3839.430.261.987543.6350.970.163.558613.5025.970.1860.6527172.0120.110.237.9568378.6019.930.1410.564934.2520.830.251.7841023.5418.300.220.9481156.9419.880.2152.43412163.8220.890.196.502合 计57.4612.4.2承压地下水资源评价承压水位于研究区中东部的含水层，受气候的影响相对潜水含水层小，所以补排量主要是人工开采、层间排泄、山区单一含水层的越流补给和迳流。研究区内迳流非常较为缓慢，迳流补给和排泄，人工开采和侧向流出是地下水主要的排泄方式，而越流补给和层间补给是承压水的主要补给来源。1995年承压水含水层组总补给量为5414.20104m3,总排泄量为59602.2104m3,均衡差为-54188.0