扩展功能
文章信息
- 李宏安, 冯雨润之, 李玲利, 贾雷
- LI Hong-an, FENGYU Run-zhi, LI Ling-li, JIA Lei
- 基于数值模拟的福州某地铁站基坑地下水控制方案研究
- Research of Foundation Pit Groundwater Control Scheme of a Metro Station in Fuzhou Based on Numerical Simulation
- 公路交通科技, 2016, Vol. 31 (2): 88-95
- Journal of Highway and Transportation Research and Denelopment, 2016, Vol. 31 (2): 88-95
- 10.3969/j.issn.1002-0268.2016.02.014
-
文章历史
- 收稿日期: 2015-01-24
2. 中国地质大学(北京) 水资源与环境学院, 北京 100013;
3. 北京地矿工程建设有限责任公司, 北京 100016
2. School of Water Resources and Environment China University of Geoscience(Beijing), Beijing 100013, China;
3. Beijing Geological Engineering Construction Co. Ltd, Beijing 100016, China
福州地区地下水埋藏较浅,且多为承压水,对地铁基坑工程的施工安全和环境安全造成极大威胁,必须采取安全稳妥的地下水控制措施。地下水控制方案设计需统筹考虑,通过解析法计算有较大的局限性,故采用数值方法对地下水位、水量的实时变化进行模拟,计算降水影响三维空间范围不同时刻水头分布。本地区存在多层淤泥质土,地层排水固结沉降量较大,控制抽水引起的沉降十分必要。
数值模拟技术在地下水控制方面应用已较为普遍,研究也较为深入。Gambolati和Freeze于1973年研究威尼斯的由多层含水层与弱透水层组成的地下水系统抽水引起的地面沉降问题时,提出两步计算模型[1]。Lewis等以此为基础提出完全耦合模型,并于1978年将其运用于威尼斯的地面沉降计算中,结果表明水头下降和地面沉降比两步计算模型较快地趋于稳定[2]。 Shearer在MODFLOW程序的基础上增加了考虑黏土层中沉降滞后于抽水的作用,将地下水流作拟三维处理,弱透水层的释水通过沉降计算后再代入含水层中,以太沙基一维固结理论为基础,开发了描述黏性土中孔隙水压力缓慢消散、滞后于含水层水头变化过程的夹层排水程序(IDP),并引入了前期固结水位的概念[3]。
杜思思等采用Processing MODFLOW软件的IBS模块对北京市平原区各含水层水位、地面沉降进行同步模拟[4]。商婷婷等采用有效应力原理的Cam-Clay模型,与GEO-SLOPE里的SEEP结合起来模拟外荷载作用下超孔隙水压力的产生和消散引起的土体固结[5]。
以美国地质调查局开发的MODFLOW商业软件为代表[6, 7, 8],太沙基有效应力原理为基础的地下水三维渗流与一维垂向固结沉降的耦合模型;骆祖江等近期提出基于比奥固结理论与土体非线性流变理论,建立地下水渗流与地面沉降三维全耦合模型[9]。金伟泽等通过上述两种方法进行对地下水渗流与地面沉降进行耦合模拟[10]。
本文以福州某地铁站为例进行模拟预测研究,首先通过抽水试验确定水文地质参数,结合场区水文地质工程条件搭建三维渗流模型,计算不同地下水控制方案下场区内水位下降情况,计算场地内外的沉降量,根据计算结果优化围护结构方案和地下水控制方案,达到控制降水引起的地层沉降的目的。
1 工程概况福州市轨道交通2号线五里亭站位于福马路与庙前街、后浦街交叉口处,呈东西向布设。现状场地地形较平坦,地面条件为中等复杂,站址周边现状为居住民房、现状道路、公共设施等,如图 1所示,布井条件有限,拆迁成本较高。降水引起的地面沉降对周边建筑物影响较大,亦可能造成部分地下管线变形甚至开裂,需根据地下水控制方案对周边环境影响的预测结果选定控制方案。
![]() |
图 1 五里亭站位置及其周边环境图 Fig. 1 Wuliting metro station and surrounding environment |
由于场区内下部承压含水层厚度较大,含水层底板埋深超过45 m,结合地下水控制要求,降水井深要求和成本控制因素,设置埋深为33 m的悬挂式隔水帷幕。车站尺寸、围护结构及工法见表 1。
车站尺寸/m | 围护结构形式 | 基坑开挖深度 | 施工工法 |
车站总长20 | 地下连续墙构筑 | 标准段基坑 | 明挖顺筑法 |
标准段宽19.7 | 悬挂式隔水帷幕 | 深度为16.2 m |
本站场区内第四系厚度约47.1 m,可分为7层(见表 2),地表以下分别为1-2杂填土,大部分地段以黏性土、碎石、块石为主回填;2-1层为灰黄色黏土,属中等压缩性土层;2-4-1层为深灰色淤泥,干强度及韧性中等,属高压缩性土层;2-4-4层为深灰色淤泥夹砂层,流塑至软塑,饱和,以黏粒为主,干强度及韧性低,属高压缩性土层;2-4-5为深灰色淤泥质细中砂,稍密至中密为主;3-1为粉质黏土,可塑至硬塑,属中等压缩性土层。3-3为(含泥)粗中砂层,饱和;3-5为深灰色淤泥夹砂,干强度及韧性低,属于高压缩性土层;3-8为浅灰色卵石层,稍密为主,饱和,间隙主要由中粗砂和砂质黏土充填。
层序 | 地层名称 | 地层厚度/m | 层底标高/m |
1-2 | 杂填土 | 1.30~4.80 | 1.45~4.84 |
2-1 | 黏土 | 2.20~7.20 | 3.43~3.80 |
2-4-1 | 淤泥 | 2.20~7.20 | -3.26~0.85 |
2-4-4 | 淤泥夹砂 | 1.20~11.20 | -23.62~-1.92 |
2-4-5 | 淤泥质细中砂 | 1.50~20.55 | -24.95~-3.62 |
3-1 | 粉质黏土 | 0.60~9.90 | -41.52~-20.60 |
3-3 | (含泥)粗中砂 | 1.10~17.30 | -42.25~-21.78 |
3-5 | 淤泥夹砂 | 1.20~2.70 | -42.85~-28.65 |
3-8 | 卵石 | 1.70~3.50 | -40.71~-33.24 |
7-1 | 花岗岩(砂土状) | 1.00~2.90 | -42.62~-42.40 |
下伏基岩层为7-1强风化花岗岩,原岩组织结构已大部分风化破坏,岩芯多呈砂土状,遇水易软化、崩解。典型地质剖面如图 2所示。
![]() |
图 2 典型地质剖面(单位:m) Fig. 2 Typical Geological section(unit:m) |
大气降雨是本区地下水的主要补给来源之一。2-4-5层淤泥质细中砂层为影响结构施工的主要承压含水层,3-3层为含泥粗中砂层,3-8层卵石承压水头高,厚度约20 m,与2-4-5层之间的相对隔水层不连续,部分地区存在“天窗”,水力联系密切。根据勘察资料,混合稳定水位埋深为0.90~3.50 m,稳定水位标高为2.86~5.75 m,地下水位变化主要受气候的控制,年变化幅度约2.0~3.0 m。地下水总体流向自西向东。
为了确定含水层的水文地质参数,查明场区内含水层之间的水力联系,于2013年6月进行了抽水试验。抽水试验场地位于福马路辅路,拟建车站的中部,共布置抽水井5眼(见图 3),沿地下水流动方向的观测井6眼,抽水试验井的结构参数见表 3。
![]() |
图 3 抽水试验孔及观测孔平面位置图 Fig. 3 Plane location of pumping test wells and observation wells |
试验分为两组多孔抽水试验:第1组对CK1-1、CK1-2和CK1-3分别进行了3个降深的抽水试验,试验目标段为12~20 m位置的淤泥质细中砂含水层,抽水试验S-lgt关系曲线见图 4(a);第2组对CK2-1,CK2-2,分别进行了2个不同降深的抽水试验,试验目标段为20~30 m位置的中粗砂含水层,抽水试验S-lgt关系曲线见图 4(b)。利用观测结果采用相应的非稳定流公式计算水文地质参数(详见表 4)。
承压非完整井,带两个观测孔渗透系数公式:

承压水条件下,带两个观测孔影响半径公式:

井号 | 含水层 | 井深/m |
井径/ mm |
过滤器 埋深/m |
过滤器 长度/m |
CK1-1 | 淤泥质细中砂 | 19.6 | 350 | 18 | 6 |
CK1-2 | 淤泥质细中砂 | 20.4 | 350 | 18 | 6 |
CK1-3 | 淤泥质细中砂 | 20.4 | 350 | 18 | 6 |
CK2-1 | 中粗砂 | 31.6 | 350 | 30 | 12 |
CK2-2 | 中粗砂 | 35.4 | 350 | 30 | 12 |
GK1-1 | 淤泥质细中砂 | 20.7 | 168 | 18 | 6 |
GK1-2 | 淤泥质细中砂 | 20.4 | 168 | 18 | 6 |
GK2-1 | 中粗砂 | 20.5 | 168 | 18 | 6 |
GK2-2 | 中粗砂 | 22.3 | 168 | 21 | 6 |
GK2-3 | 中粗砂 | 20.5 | 168 | 18 | 6 |
GK2-4 | 中粗砂 | 25.3 | 168 | 24 | 12 |
注:第2组井(CK2-1,CK2-2,GK2-1,GK2-2,GK2-3,GK2-4)针对上部淤泥质细中砂层含水层进行了分层止水处理。 |
![]() |
图 4 抽水试验S-lgt关系曲线 Fig. 4 Curve of relationship between S and lgt(pumping tests) |
渗透系数K/(m·d-1) | 影响半径R/m | ||
第1组 | 1 | 11.64 | 138.66 |
2 | 7.04 | 50.59 | |
第2组 | 3 | 26.14 | 58.23 |
4 | — | — | |
5 | — | — | |
6 | 34.05 | — | |
7 | 22.85 | 212.67 | |
8 | 31.28 | 228.46 | |
综合 建议值 | 第1组 | 7.04 | 80~120 |
第2组 | 22.85 | 200~240 |
式中,S1为抽水稳定后1号观测孔水位降深;S2为抽水稳定后2号观测孔水位降深;r1为1号观测孔到主孔的距离;r2为2号观测孔到主孔的距离; h1,h2为抽水稳定时1号观测孔和2号观测孔距含水层底板的水位高程;Q为抽水孔稳定流量。为了解抽取地下水地面沉降情况,本次抽水试验对可能影响到的区域范围内布置地面沉降观测点,结合场地实际情况共布置22个地面沉降观测点(沿着抽水试验影响范围均匀布置),在抽水试验期间定期进行沉降观测。沉降观测时间为:7月16日—7月25日。
周边10 m范围内监测点平均累计沉降18.78 mm,周边10~30 m范围内监测点平均累计沉降13.19 mm,周边30~50 m范围内监测点平均累计沉降11.59 mm,周边50 m处监测点累计沉降0.46 mm。仅仅短暂的单井抽水试验即引起如此范围和量级的地面沉降,如实施大规模、长时间的车站基坑降水,而不采取其他控制措施,对周边环境的破坏影响可想而知。
3 地下水渗流模型的建立地下水渗流模型可预测各种条件影响下(比如时间变化,边界条件变化,不同部位的源、汇变化等)的渗流场变化情况,这些变化条件与实际工程施工相对应,可针对性较强地增减泵量、调整开泵部位和区域;以此为基础进一步研究地下水位以下各地层由于水位降低有效应力增加所产生的地面沉降变化规律。
3.1 模拟软件介绍本次模拟选用Visual MODFLOW数值分析软件。Visual MODFLOW是目前国际上最为流行的三维地下水流和溶质运移模拟的标准可视化专业软件。由于其程序结构的模块化、有限差分离散的简单化和求解方法的多样化等优点,被广泛用来模拟井流、河流、排泄、蒸发和补给对非均质和复杂边界条件的水流系统的影响。
3.2 数学模型多孔介质承压含水层三维非稳定流数学模型方程如下所示[11]:

式中,Ω为渗流区域;h为地下水系统的水位标高;K为含水介质的水平渗透系数;Kz为含水介质垂向渗透系数;ε为含水层的源汇项;h0为系统的初始水位分布;S为自由面以下含水层储水率;Γ0为渗流区域的上边界,即地下水的自由表面;μ为潜水含水层在潜水面上的重力给水度;p为潜水面的蒸发和降水入渗强度等;Γ1为已知水位边界;h1为已知边界水位值;Γ2为渗流区域的流量边界;Kn为边界面法线方向的渗透系数,n为边界面的法线方向;q为Γ2边界的单位面积上的流量,流入为正,流出为负,隔水边界为0。
3.3 模型范围和边界条件场地距离自然边界较远,为减少地下水边界效应的影响,由基坑中心向四周扩展2倍的影响范围(500 m),即模拟区为1 000 m×1 000 m的矩形区域。在模型边界处设定定水头边界,并在模型的东、西边界上设置一定的水头差以模拟天然状态下地下水的渗流。
结合勘察资料和该地区水文地质条件,将模拟区的第四系地层概化处理为5层,地层概化参数见表 5。
层序 | 含水层性质 | 主要岩性 |
底板标 高/m |
渗透系数/ (m·d-1) | 孔隙比 |
1 | 隔水层 | 杂填土、黏土、淤泥 | -5.6 | 0.5 | 1.5 |
2 | 承压含水层 | 淤泥质细中砂 | -17.1 | 8.5 | 1.6 |
3 | 隔水层 | 粉质黏土、淤泥夹砂 | -20 | 0.05 | 0.75 |
4 | 承压含水层 | (含泥)粗中砂 | -34.5 | 25 | 0.6 |
5 | 承压含水层 | 卵石 | -40.15 | 40 | 0.5 |
6 | 隔水帷幕 | 钢筋混凝土 | — | 1e-05 | — |
模型采用六面体网格剖分,在水平方向上采用非等距矩形网格剖面(基坑开挖区域附近网格加密)。综合考虑隔水帷幕的设置、边界条件、水文地质条件、抽水井和观测井的布置,平面上剖分为60行、107列,加密区最小单元格的面积为2×2 m2,非加密区域单元格面积约为20×20 m2;垂向上共划分为5个参数分区。剖分网格见图 5。
![]() |
图 5 三维离散网格剖分示意图(单位:m) Fig. 5 Three-dimension subdivision of discrete model nets(unit:m) |
根据勘察报告中的地下水水位,按照内插法和外推法获得含水层的初始水位,考虑整体流场赋虚拟水位值,基坑降水模型源汇项比较简单,水面可近似作为平面考虑,给整个模型赋一个水头值,然后通过稳定流运行得到初始流场。降水工程中常需要疏干某一含水层,而MODFLOW的计算原理决定了它无法良好地解释含水层疏干条件下的渗流情况,因此在确定模拟期时,需要参考相关工程经验,将总模拟期控制在水位下降至疏干的时间段内。本次研究模拟抽水30 d内的变化,将30 d分为10个应力期,时间步长为3。
4 方案研究及预测分析 4.1 方案比选本区域地层软弱,以淤泥、黏土、淤泥质砂为主。随着降水作业的进行,有效应力增大导致的土体压缩将会非常明显,因此本次地下水控制设计采用悬挂式帷幕,并将降水井限定在帷幕范围内,这将大大降低对周边环境的影响。在这一前提下降水引起的沉降量采用了《建筑基坑支护技术规程JGJ120—2012》中相关公式结合数值模拟结果进行计算。
固结沉降的计算方法为分层总和法,公式为:

式中,s为降水引起的地层变形量;ψw为沉降计算经验系数,应根据地区工程经验取值,无经验时,宜取ψw=1;Δhi为第i层土的厚度;Δσ′zi为降水引起的地面下第i土层中点处的附加有效应力,对黏性土,应取降水结束时土的固结度下的附加有效应力;Esi为第i层土的压缩模量,应取土的自重应力至自重应力与附加有效应力之和的压力段压缩模量值。
帷幕内抽水降压采用多点小流量分散降压抽水设置降水井,沿车站主体布置两排共42口降压井,在车站东西端头和中心各设置一口观测井,见图 6,模拟抽水30 d水位的变化情况。降水井间距 9 m,井深30 m,模型内设置为分层降压,其中20口井只对4层降压,22口对2层和4层降压。模型中抽水井出水量根据各方案实际情况有所调整,见表 6。
![]() |
图 6 降水井平面布置图 Fig. 6 Plane layout of dewatering wells |
2层单井出水量/(m3·d-1) | 4层单井出水/(m3·d-1) | ||
1~15 d | 15~30 d | 1~30 d | |
方案1 | 10 | 5 | 125 |
方案2 | 10 | 5 | 112 |
方案3 | 20 | 20 | 64 |
本方案围绕车站主体设置悬挂式帷幕,井深30 m,隔水帷幕深度为33 m,基坑不进行封底处理,该方案垂向水头分布见图 7。由于3层隔水层车站右线部分地段厚度较小,甚至出现天窗,出于安全考虑,2层和4层均降至基坑结构底板以下0.5 m,即要求降深为13.78 m。抽水30 d后,时间和降深的关系见图 8。
![]() |
图 7 方案1垂向水头分布图(单位:m) Fig. 7 Vertical distribution of hydraulic heads of scheme 1(unit:m) |
![]() |
图 8 方案1时间与降深关系曲线 Fig. 8 Curves of relationship between drawdown and time of scheme 1 |
帷幕范围内抽水造成的沉降约270~310 mm,最大沉降量位于车站中段,约330 mm;帷幕外地面最大沉降量为140 mm(图 9)。
![]() |
图 9 方案1地面沉降等值线(单位:mm) Fig. 9 Land subsidence contour lines of scheme 1(unit:mm) |
本方案围绕车站主体设置悬挂式帷幕,井深30 m,隔水帷幕深度为38 m,基坑不进行封底处理,该方案垂向水头分布见图 10。含水层2层和含水层4层的要求降深均为13.78 m(同方案1),即车站主体结构底板以下0.5 m。抽水30 d后,时间和降深的关系见图 11。
![]() |
图 10 方案2垂向水头分布图(单位:m) Fig. 10 Vertical distribution of hydraulic heads of scheme 2(unit:m) |
![]() |
图 11 方案2时间与降深关系曲线 Fig. 11 Curves of relationship between drawdown and time of scheme 2 |
帷幕范围内抽水造成的沉降约270~310 mm,最大沉降量位于车站中段,约320 mm;帷幕外地面最大沉降量为130 mm(图 12)。
![]() |
图 12 方案2地面沉降等值线(单位:mm) Fig. 12 Land subsidence contour lines of scheme 2(unit:mm) |
本方案围绕车站主体设置悬挂式帷幕,井深30 m,隔水帷幕深度为33 m,采用三轴搅拌桩在基坑结构底板以下做旋喷封底处理,封底厚度3 m,该方案垂向水头分布见图 13。含水层4层水位降至满足抗突涌要求(降深7.48 m)。含水层2层封底以上部分疏干,封底以下部分水位降至满足抗突涌要求(降深11.28 m)。抽水30 d后,时间和降深的关系见图 14。
![]() |
图 13 方案3垂向水头分布图(单位:m) Fig. 13 Vertical distribution of hydraulic heads lines of scheme 3(unit:m) |
![]() |
图 14 方案3时间与降深关系曲线 Fig. 14 Curves of relationship between drawdown and time of scheme 3 |
帷幕范围内抽水造成的沉降约170~210 mm,最大沉降量位于车站中段,约220 mm;帷幕外地面最大沉降量为80 mm(图 15)。
![]() |
图 15 方案3地面沉降等值线(单位:mm) Fig. 15 Land subsidence contour lines of scheme 3(unit:mm) |
3个比选方案的各项参数统计见表 7。
方案 |
帷幕深 度/m |
旋喷 封底 |
坑内最大 沉降量/m |
坑外最大 沉降量/m | 基坑总涌水量/(m3·d-1) | |
1~15 d | 15~30 d | |||||
方案1 | 33 | 无 | 320 | 140 | 5 470 | 5 360 |
方案2 | 38 | 无 | 310 | 130 | 4 924 | 4 814 |
方案3 | 33 | 有 | 220 | 80 | 3 128 | 3 128 |
通过3个方案的比较,可以发现方案2的悬挂式帷幕由方案1的33 m加深至38 m,改变了坑外水的补给渗流路径,导致涌水量减少,减少了抽水设备的工作强度,但是对于基坑周边地面沉降的控制效果无明显增强;方案3通过旋喷封底形成了一个等效隔水层,2层细中砂在封底内的部分厚度仅为4 m,极易疏干;同时由于旋喷封底水泥土重度略大,4层含泥粗中砂层的降压要求也有所降低,使得降排水量减少,降水引起的地面沉降量也大为减小,最大沉降量为80 mm,比方案2少了50 mm,效果显著。
5 结论(1)对比帷幕内外地下水位变化引起的沉降结果表明,设置隔水帷幕对沉降控制的效果明显。但是悬挂式帷幕在深度已超过降水井井深的情况下,仅仅加大帷幕深度,对周边沉降的控制和对基坑涌水量的减小效果不明显,远远不如落底式隔水帷幕。
(2)方案3对基坑底部做了旋喷封底处理,承压水压力水头降低减小,大大减少了基坑涌水量,对周边环境的影响相对其他两个方案也最小,因此建议采用方案3。
(3)方案3在帷幕外侧产生80 mm的最大沉降,周边环境(特别是地下管线)也难以承受,因而需采取地下水回灌措施以减小地铁基坑外侧承压水头降深值,把降水引起的地面沉降量控制在环境安全允许范围内。
(4)建议在本次研究的基础上建立降水-回灌的数值模拟模型,研究在坑外进行地下水回灌条件下降水对周边环境的影响,并确定可靠的降水-回灌方案。
(5)通过必要的补充勘察,进一步查明两个承压含水层之间“天窗”范围,研究存在“天窗”水力联系的降水-回灌的数值模拟模型,使数值模拟的分析研究成果更加切合工程实际。
[1] | GAMBOLATI R,FREEZE A.Mathematical Simulation of the Subsidence of Venice:Theory[J]. |
[2] | LEWIS R W,SCHREFLER B.A Fully Coupled Consolidation Model of the Subsidence of Venice[J]. |
[3] | SHEARER T R.A Numerical Model to Calculate Land Subsidence,Applied at Hangu in China[J]. |
[4] | 杜思思.海河平原地下水与地面沉降模型模拟研究[D].北京:中国地质大学,2011.DU Si-si.Study on Simulation of Model of Ground Water and Land Subsidence in Haihe River Plain[D].Beijing:China University of Geosciences,2011. |
[5] | 尚婷婷.天津市塘沽区城市建设引发地面沉降数值模拟研究[D].北京:中国地质大学,2007.SHANG Ting-ting.Research on Simulation of Land Subsidence Caused by Urban Construction in Tanggu District of Tianjin[D].Beijing:China University of Geosciences,2007. |
[6] | 阚京梁,罗立红.Processing Modflow模型在预测地面沉降中的应用[J].铁道工程学报,2010,27(2):27-31.KAN Jing-liang,LUO Li-hong.Application of Processing Modflow Model in Prediction of Ground Subsidence[J].Journal Railway Engineering Society,2010,27(2):27-31. |
[7] | 李英,何钟泽,严桂华,等.武汉二元结构地层基坑降水及其地面沉降研究[J].岩土工程学报,2012,34(增1):767-772.LI Ying,HE Zhong-ze,YAN Gui-hua,et al.Excavation Dewatering and Ground Subsidence in Dual Structural Stratum of Wuhan[J].Chinese Journal of Geotechnical Engineering,2012,34(S1):767-772. |
[8] | 周念清,唐益群,娄荣祥,等.徐家汇地铁站深基坑降水数值模拟与沉降控制[J].岩土工程学报,2011,33(12):1950-1956.ZHOU Nian-qing,TANG Yi-qun,LOU Rong-xiang,et al.Numerical Simulation of Deep Foundation Pit Dewatering and Land Subsidence Control of Xujiahui Metro Station[J].Chinese Journal of Geotechnical Engineering,2011,33(12):1950-1956. |
[9] | 骆祖江,曾峰,李颖.地下水开采与地面沉降控制三维全耦合模型研究[J].吉林大学学报:地球科学版,2009,9(6):1080-1086.LUO Zu-jiang,ZENG Feng,LI Ying.Study on Study on Three-dimensional Full Coupling Model of Groundwater Exploitation and Land-subsidence Control[J].Journal of Jilin University:Earth Science Edition,2009,9(6):1080-1086. |
[10] | 金玮泽,骆祖江,陈兴贤,等.地下水渗流与地面沉降耦合模拟[J].中国地质大学学报:地球科学版,2014,39(5):611-619.JIN Wei-ze,LUO Zu-jiang,CHEN Xing-xian,et al.Coupling Simulation of Groundwater Seepage and Land Subsidence[J].Journal of China University of Geosciences:Earth Science Edition,2014,39(5):611-619. |
[11] | 薛禹群.地下水动力学[M].北京:地质出版社,1986.XUE Yu-qun.Groundwater Dynamics[M].Beijing:Geological Publishing House,1986. |