含水层水文地质

含水层水文地质（精选8篇）

含水层水文地质 第1篇

运裕煤业有限责任公司设计生产能力90Mt/a, 地处太行山西北部, 属丘陵地貌, 地表多为黄土覆盖, 井田内发育有两条季节性河流及一座水库, 河流雨季流量较大, 大雨时常见有山洪暴发, 水位及流量观测数据缺乏。洪水川位于井田北界, 河水由西向东径流。南峪河位于井田东部, 横穿井田东南部由南向北流入洪水川后往北注入松溪河。松溪河上游的郭庄水库为井田的东南边界。河流附近区域15#煤覆盖层厚度仅为100-140m, 如何解决地表河流及冲积含水层对矿井生产安全造成的不利影响是运裕煤矿面临的巨大挑战。

2 松散层探查

107工作面位于首采区东部, 直接受到南裕河影响, 且工作面南部地表发育有沙沟, 地势及地表高程相对较低;而根据煤层底板等值线图, 工作面南部煤层底板标高相对较高。所以, 该面南部处于松散含水层的汇水区域, 且覆岩厚度相对较薄, 煤层回采后顶板导水裂缝带将波及松散含水层。故松散含水层的地层结构及水文地质条件对该面开采工艺及防排水能力配备均具有重要影响。

为了探明松散含水层地层结构及水文地质条件, 该矿在沙沟处施工了松散层探查孔, 开孔孔径Φ168mm, 松散层段下入Φ146mm滤管, 终孔孔径Φ127mm, 终孔深度为30.1m, 进入基岩15.5m。探明该处松散层厚度为14.6m, 风化带厚度为8.7m, 松散层顶部为1.2m厚的黄土层, 其下为卵石层、砂土层、粘土层互层结构。其中, 粘土层自上而下发育有3层, 厚度分别为2.8m、0.8m、0.48m, 总计4.08m。这三层粘土层虽然厚度均较小, 但是其较强的遇水膨胀性, 可以起到快速闭合采后裂隙、阻隔松散层进入采掘空间的作用, 故对107面的防治水工作具有较大意义。

3 水文地质参数计算

松散层探查孔进行了单孔稳定流抽水试验, 试验结束后, 采用潜水完整井裘布依公式[1]计算松散含水层水文地质参数:

式中K含水层渗透系数 (m/d) ;

Q抽水井流量 (m3/d) ;sw抽水井中水位降深 (m) ;R影响半径 (m) ;

H潜水含水层厚度 (m) ;

h潜水含水层抽水后的厚度 (m) ;rw抽水井半径 (m) 。

通过探查孔探明该处松散含水层厚度1.67m, 水位标高+888.683m, 单位涌水量0.5801l/s.m, 渗透系数41.5m/d, 影响半径21.4m, 为富水性中等含水层, 连通性较好。

4 矿井涌水量预测

根据107面水文地质条件分析可知, 该工作面充水水源主要有三类, 即南裕河、松散含水层、太原组灰岩含水层。由于采用限制采厚综采开采工艺, 南裕河水不会直接进入采掘空间, 仅以渗流方式补给松散含水层;太原组灰岩含水层富水性弱, 且根据周边已经回采工作面实测资料, 采用比拟法可得107面回采时, 太原组灰岩含水层贡献的正常涌水量为40m3/h, 最大涌水量为70m3/h;松散含水层贡献的矿井涌水量则需要专门计算, 本次采用“大井法”[2]进行预

4.1 基本原理

矿井排水时, 其周围会形成以巷道系统为中心的具有一定形状的降落漏斗, 这与钻孔抽水时在钻孔周围形成降落漏斗的情况相类似, 因而可以将巷道系统分布的范围假设为一个理想的大井。这个假设大井的圆形断面积, 与巷道系统分布的面积相当。因此, 可以直接利用地下水动力学的公式来计算巷道系统的涌水量。

4.2 水文地质概念模型

4.2.1 充水岩层的边界性质和形态

周边界的概化:本次计算区域呈矩形, 长度a=430m, 宽度b=150m, 周边界为补给边界。

内边界概化及引用半径r0的确定:内边界是指井、钻孔或巷道系统的进水边界。对于15107工作面来讲, 内边界是矩形, 引用半径rO适用公式为rÁ4ab, 又因ab0.35, 故1.14, 则Ár 165m。

4.2.2 充水岩层的介质性质

根据探查孔资料, 松散层厚度为1.67m, 将其概化为理想的均值、等厚含水层。

4.3 计算公式

松散含水层属于潜水含水层, 故选用大井法潜水公式进行涌水量计算:

式中:Q预计的矿井涌水量, m3/d;

K含水层渗透系数, m/d;

H潜水含水层厚度, m;

S由于矿井排水而产生的水位降低值, m;

R含水层抽水时得出的影响半径, m;

r0假想大井的半径 (或称引用半径) , m;

R0矿井排水的引用影响半径, R0=R+r0, m。

通过计算可得, 107面回采时, 松散层含水层水进入矿井的涌水量为125 m3/h。

故107面回采过程中, 因为松散含水层及太原组灰岩含水层组合形成的矿井正常涌水量为165 m3/h, 最大涌水量为195 m3/h。根据此计算结果, 安排该工作面配备不小于400 m3/h的排水能力。

5 结论建议

本文介绍了运裕煤矿107工作面南部沙沟松散层含水层水文地质条件探查方法, 探明了松散层岩性组合情况, 计算了松散含水层水文地质参数。在上述工作基础上, 运用大井法计算了松散层含水层产生的矿井涌水量, 运用比拟法计算了太原组灰岩含水层产生的矿井涌水量。在该面实际生产过程中, 稳定涌水量约为150 m3/h, 最大涌水量约为200 m3/h, 与计算结果基本一致, 且由于提前配备了不小于400 m3/h排水能力的排水设备, 保证了该工作面安全回采。

摘要：阐述了松散层地质、水文地质条件探查在水体下采煤工作中的重要性, 以运裕煤矿107工作面南部沙沟松散含水层水文地质条件探查工作为基础, 总结了松散含水层水文地质条件探查技术方法、工作流程, 提供了涌水量计算方法, 为水体下采煤安全性评价提供了技术支持。

关键词：松散含水层,探查,涌水量计算

参考文献

[1]虎维岳.矿山水害防治理论与方法[M].北京:煤炭工业出版社, 2005:65-78.

邢东矿主要含水层水质特征分析 第2篇

邢东矿主要含水层水质特征分析

文章介绍了邢东井田概况,并重点对邢东矿2001～67份不同含水层的水质资料通过绘制水质图进行了系统分析总结,得出了各主要含水层的水质特征,从而为矿井快速判别出水水源积累了经验.

作 者：高春芳 GAO Chun-fang  作者单位：冀中能源金牛股份公司邢东矿,河北邢台,054001 刊 名：河北煤炭 英文刊名：HEBEI COAL 年，卷(期)：2009 “”(1) 分类号：P641.3 关键词：含水层   水质   特征   分析  

含水层水文地质 第3篇

罗园井田地处淮河中游, 属淮河冲积平原, 区内地形平坦, 西高东低, 地面标高一般在+24~+28m。勘查区东侧有西淝河流入淮河。

本区属过渡带气候, 为季风温暖带半湿润气候, 季节性明显, 夏季炎热, 冬季寒冷。

多年平均降雨量为902.9mm, 最大1 313.5mm, 最小389.8mm, 日最大降雨量173.1mm, 小时最大降雨量75.3mm, 降雨多集中在6、7、8三个月, 约占全年降雨量的40%;每年初雪一般在11月上旬, 终雪在次年3月中旬, 雪期72~127天, 年平均蒸发量1 695.5mm (水面) , 蒸发量大于降水量, 潮湿系数近似0.5, 最大相对湿度78%, 最小相对湿度10.14%, 平均74%。

2 区域构造

淮南煤田位于秦岭纬向构造带南亚带的北缘, 东与新华夏郯城~庐江断裂反接, 西连周口凹陷, 北接蚌埠隆起, 南邻合肥凹陷。

淮南煤田为一复向斜构造, 由于南北应力的挤压作用, 促使淮南复向斜主体构造形迹呈东西向展布, 并在复向斜南北两翼发育了一系列压扭性的逆冲断层, 造成了两翼的迭瓦式构造, 使南翼部分地层直立倒转, 由于新华厦系构造的干扰, 使主体构造形迹稍有变动, 略呈北西西~南东东向展布, 这一纬向构造带经过燕山运动的改造, 形成了现在的构造轮廓, 不仅表现褶皱、断裂, 并伴有岩浆活动。

3 物性特征

3.1 第一含水层 (组)

第一含水层 (组) 埋深一般在5.0m~25m左右, 砂层发育较好, 主要由细砂、粉砂、粘土质砂组成, 从曲线的形态和幅值可分为二种类型:

1) 从曲线的形态分一含上、下部为厚层砂, 视电阻率曲线幅值反映较高, 一般在30Ω.M左右, 上部砂层厚度为6m左右, 下部砂层厚度为8m左右。其物性反映特征见图1;

2) 一含砂层发育较好, 可分为上、中、下三段, 视电阻率值一般在30Ω.M左右。下部砂层一般厚度为5m左右, 质纯, 粒度均一, 曲线反映圆滑。

3.2 第二含水层 (组)

第二含水层 (组) 埋深一般在45m~65m左右, 砂层发育较稳定, 一般有一厚层砂或上、下部二层砂组成, 曲线异常反映明显, 界面清晰, 主要以细砂、粉砂、粘土质砂组成, 从曲线的形态和幅值可分为2种类型:

1) 二含砂层发育好, 有一厚层砂组成, 视电阻率曲线幅值一般在30Ω.M~40Ω.M左右, 含水层厚度在15m左右;

2) 二含上、下部砂层发育较好, 视电阻率曲线幅值一般在30Ω.M左右, 中部多为厚层粘土。

3.3 第三含水层 (组)

第三含水层埋深一般在80m~125m左右, 砂层发育较好, 一般有三层或四层厚砂组成, 曲线异常明显, 界面清晰, 岩性主要以细砂和粉砂为主, 夹有少量中砂, 从曲线的形态和幅值可分为4种类型:

1) 三含砂层发育好, 有四厚层砂组成, 视电阻率曲线幅值一般在30Ω.M左右, 含水层厚度在40m左右;

2) 三含上、中、下部砂层发育较好, 视电阻率曲线幅值一般在50Ω.M左右, 中部多为厚层粘土, 含水层厚度在25m左右;

3) 三含上、中部砂层发育较好, 视电阻率曲线幅值一般在40Ω.M左右, 含水层厚度在25m左右;

4) 三含下部砂层发育较好, 上部为粘土质砂、粘土互层, 视电阻率曲线幅值一般在30ΩM左右。

4 含水层水文地质条件特征

4.1 第一含水层 (组)

一含地层时代属第四系全新统, 为冲积平原河漫滩相~泛滥相沉积, 含水层类型为孔隙潜水~弱承压水含水层, 砂层顶界面埋深2.10m~14.80m, 一般3m~6m, 底界面埋深22.90m~30.65m, 一般24m~27m, 岩性主要为黄色、青灰色及浅灰绿色粉砂及粘土质砂, 局部发育薄层细砂, 其中夹粘土0~6层, 砂层厚度为5.40m~18.50m, 一般12m~15m, 平均厚度为12.34m, 在勘探区内砂层西北部较厚, 东南部较薄。

4.2 第二含水层 (组)

二含地层时代属第四系更新统中下部, 为河床相~河漫滩相沉积, 含水层类型为孔隙承压含水层, 顶界面埋深41.60m~50.80m, 一般44m~48m, 底界面埋深59.85m~71.60m, 中部埋藏较浅, 向两边逐渐变深。岩性主要为黄色、灰黄色、灰绿色中细砂及粘土质砂, 砂层主要成分为石英、长石, 含少量暗色矿物, 砂层厚度10.35m~25.20m, 一般厚15m~19m, 平均厚度为16.67m。砂层厚度较薄, 厚度变化不明显, 但总体趋势是勘探区东部和西部偏厚, 中部变薄。

4.3 第三含水层 (组)

本层为河漫滩相~湖泊相沉积的孔隙承压含水层, 含水砂层顶界面埋深74.00m~83.40m, 底界面埋深117.00m~138.90m, 岩性主要为灰黄色、灰绿色中细砂为主, 次为粘土质砂, 主要成分为石英、长石, 含少量暗色矿物。砂层厚度23.15m~59.80m, 一般厚27m~38m, 平均厚度为35.44m, 在勘探区内砂层东北和西南部较厚, 东南部、中部及西北部较薄。

参考文献

[1]李学富.桃园煤矿松散底部含水层水文地质特征分析.安徽科技, 2010 (7) .

[2]王红喜.鹤壁六矿主要含水层水文地质特征及水质分析.中州煤炭, 2004 (1) .

[3]高明志.山东沂源下寒武统馒头组含水层水文地质特征.山东地质, 1997 (2) .

含水层水文地质 第4篇

河渠控制下潜水含水层中的蒸发规律与计算

在潜水蒸发和河渠水平渗透共同作用下,对潜水非稳定渗流模型,利用Boussinesq第一线性化方法,通过Laplace变换,导出了该模型的解析解;对时变的潜水蒸发强度,利用卷积原理,给出计算潜水水位变动的.解析公式;建立了计算潜水蒸发强度的递推公式;根据模型解,研究了潜水蒸发引起的潜水位变化过程与变化规律.认为在河渠边界控制下的潜水含水层中,蒸发作用将形成具时间滞后特征的潜水位变化过程,河渠边界对潜水位变化过程的影响幅度具有非线性时变特征;如果忽略边界条件所进行的评价,其计算误差将随时间的延长而迅速增大.以安徽淮北平原某河流-潜水含水层系统为例,阐述上述方法过程与步骤.

作 者：陶月赞 作者单位：合肥工业大学土建学院,合肥,230009刊 名：自然科学进展 ISTIC PKU英文刊名：PROGRESS IN NATURAL SCIENCE年，卷(期)：16(5)分类号：N1关键词：潜水非稳定流 潜水蒸发 递推公式 非线性时变 滞后

含水层水文地质 第5篇

关键词：含水层,水文地质特征,朱仙庄矿

我国煤矿水文地质条件是极为复杂的, 地表及以下往往是由多含水层构成的复杂系统。一方面, 断层、陷落柱以及含水层间的越流补给形成了地下水系统复杂的水循环模式。另外一方面, 煤层开采造成地下岩层应力的重新分布和岩层的破裂损伤, 这种破裂损伤极大地改变了围岩的渗透性, 使地下主要充水含水层受到不同程度的扰动, 导致生产过程中出现涌水甚至突水, 从而改变了含水层原来的地下水运行状态, 破坏了原有的地下水平衡。水害是矿井的一大灾害, 地下水的防治是煤矿安全生产的一项重大内容, 它也是影响煤矿生产的一个重要因素。煤炭资源的开采受水害威胁严重, 尤其是随着开采深度、开采强度的、开采规模的增大和扩大。我国华北型石炭二叠系煤田地层多为奥陶纪石灰岩, 岩溶发育, 且又往往与断层、裂隙等相连通, 富含承压水, 水害事故频繁发生, 严重的威胁煤矿的安全生产, 并已造成巨大的经济损失。研究矿井涌水对煤矿生产的危害, 一般从矿井水源和涌水通道两方面着手进行处理。治理矿井涌水的关键问题是快速、正确判定涌水水源的形成。对含水层水文地球化学特征的研究能有助于分析含水层的水力联系和判定突水水源。

一、井田位置及地质概况

朱仙庄煤矿位于宿州市东南13km处, 属宿州市埇桥区管辖, 西北距淮北市64km。其边界:南与芦岭煤矿相邻, 北、东、西均以煤层露头为界, 深部至各煤层-800m水平线;南北走向长9km, 东西倾向宽1.5~5.8km, 面积21.555km2。井田位于宿东向斜的北段。宿东向斜位于徐宿弧形构造南端, 西寺坡逆断层和东三铺逆断层所夹持的断块之间。轴部为二叠系煤系地层, 四周为石炭系、奥陶系灰岩所包围, 东翼因东三铺逆断层切割而使寒武、奥陶系灰岩压覆到煤系地层之上。井田均为新生界松散层所覆盖, 经钻孔揭露地层有奥陶系、石炭系、二叠系、侏罗系和第四系。

二、主要充水含水层水文地质特征分析

矿区主要含水层有第三、四系含水层 (组) ;侏罗系第五含水层 (组) ;二叠系含水层 (段) ;太原组灰岩岩溶含水层 (段) (八含) ;奥陶系灰岩岩溶含水层 (段) (九含) , 其中第三、四系含水层 (组) 按其岩性及含水性自上而下分为四个含水层 (分别为第一含水层、第二含水层、第三含水层、第四含水层) 。

(一) 第三、四系含水层 (组)

第三、四系为新生界松散层, 总厚218~261m, 按其岩性及含水性自上而下分为四个含水层 (组) 。第一含水层 (组) :由浅黄色粉砂及粘土质砂夹薄层砂质粘土组成, 水质类型为HCO3NaMg型, 该层水可作为供水水源。第二含水层 (组) :由3~6层灰黄色粉砂及细砂组成, 水质HCO3NaMg型, 是本区主要供水水源。第三含水层 (组) :由灰黄色、灰绿色含泥质的粉、细、中砂组成, 水质类型为SO4ClNaMg型。第四含水层 (组) :本组含水层岩性、厚度变化受古地形控制, 而富水性大小又与岩性密切相关, 在区内四含为一封闭的盆地沉积。第四系松散层一、二、三含水层 (组) , 由于受第三系第三隔水层的阻隔, 对矿井充水无直接影响, 而第四含水层组与煤系地层直接接触, 是矿井的主要补给水源之一。

(二) 侏罗系第五含水层 (组)

本组分布在朱仙庄矿东北部, 经物探和钻探查明井田以北有的塔桥断层, 使下盘五含成为孤立块段, 五含的影响面积由原来的500km2缩小为9km2, 其中矿内有2.8km2。五含与四含及煤系地层均为不整合接触, 属山麓洪积相沉积, 砾石成份以灰岩较多, 主要为灰岩碎块, 分选差, 胶结物为紫红色泥岩, 厚度受古地形控制, 浅部和西部大, 向深部和东部逐渐变薄和尖灭, 剥蚀面大致与8煤层平行, 倾向东北。五含岩溶极为发育, 但不均匀, 充填有泥质、方解石和石膏等物。五含的富水性决定于岩溶的发育程度, 由于岩溶发育很不均匀, 因此不同区段的富水性强弱不等, 富水性随着砾岩厚度由南向北增加, 富水性增强。

(三) 二叠系含水层 (段)

煤系由泥岩、粉砂岩及砂岩组成, 裂隙不发育, 泥浆消耗量小, 含水性弱, 含水层不易划分, 现依据煤层顶、底板岩性将二叠系地层划分为二个含水层 (组) 。下石盒子组第六含水层 (段) :位于二叠系上石盒子组上部, 界于四煤与三煤间, 含水层为细~中粒砂岩, 平均厚90m, 俗称K3砂岩含水层 (组) 。本含水层富水性弱, 距下部七、八开采煤层较远, 对开采影响不大。下二叠系第七含水层 (段) :位于二叠系下部七、八、十煤层附近, 岩性为中、细粒砂岩, 平均厚70m左右。砂岩水的富水性取决于砂岩裂隙的发育程度, 矿井南北略有差异。由于松散层第四含水层直接覆盖在该含水层露头, 当采动破坏, 砂岩水直接进入矿井, 同时四含水通过砂岩裂隙顺层补给进入矿井, 所以当井下工作面出水后, 煤层顶底板砂岩水首先下降。同时四含水位也相应下降, 说明煤层顶底板砂岩裂隙水, 不仅是矿井充水的直接水源, 同时也是四含水进入矿井的主要补给通道。

(四) 太原组灰岩岩溶含水层 (段) (八含)

该含水层总厚140m左右, 含灰岩11~12层, 灰岩总厚62m左右, 约占该组的44%, 各层灰岩间有一定的泥岩、粉砂岩隔水层。灰岩的富水性取决于岩溶裂隙发育程度, 根据区内四个孔的抽水资料, 一般浅部露头区富水性较强;深部-300m以下较弱, 有明显的垂直分带规律。水质为HCO3SO4NaCa型, 矿化度1.04~1.79g/l, 富水性中等。该含水层距10煤层60m左右, 由于该含水层水位较高, 部分受断层影响, 隔水层变薄地带, 对10煤层的开采都有较大的威胁。

(五) 奥陶系灰岩岩溶含水层 (段) (九含)

该层区域探明厚度大于500m, 分布在本矿两翼, 远离煤层, 对矿井开采无直接充水影响, 据抽水资料, 水质为SO4ClNaCa或CaMg型。富水性与岩溶裂隙发育程度密切相关, 但总的来说, 该含水层富水性不均一, 局部富水性强。

三、充水因素分析

(一) 矿井充水水源

本矿地层中有多个不同含水层 (组、段) , 但它们在垂直剖面上相应也有多个不同的隔水性良好的隔水层 (组、段) 所阻隔, 不同 (组、段) 的地下水对矿井的充水影响有明显的不同。地表水和新生界松散层第一、二、三含水层 (组) 地下水由于受第三隔水层 (组) 所阻隔, 它们一般对矿井充水无影响。其主要充水水源有如下几个:1) 新生界松散层:第四含水层 (组) 矿内普遍分布, 厚度0~49m。岩性复杂, 其富水性弱~中等。侏罗系第五含水层在本矿东北部发育, 富水性较强。四含、五含地下水可沿浅部基岩风化带裂隙、煤层采空冒裂带裂隙和顺煤层进入矿井, 在留有防水煤柱的情况下是浅部煤组开采的补给水源。2) 各主采煤层顶、底板砂岩裂隙含水层:是矿井充水的直接含水层, 由于砂岩裂隙发育具有不均一性, 一般富水性较弱。地下水处于封闭~半封闭环境, 补给条件差。3) 煤系内断层:在天然状态下断层中泥质较多, 断层影响带富水性弱, 导水性差。当断层落差较大时, 形成太原组石灰岩 (或奥灰) 岩溶裂隙水与主采煤层对口或间距缩小, 采煤时有发生突水的可能性。4) 灰岩岩溶裂隙水:在正常情况下下部灰岩岩溶裂隙水对可采煤层开采无直接充水影响, 但本矿逆断层较发育, 且断距较大, 使可采煤层与灰岩岩溶裂隙多处对口接触。当井巷工程遇到导水断层或陷落柱时, 灰岩岩溶裂隙水可能对矿坑产生直接充水或发生“底鼓”。由于灰岩水水压大, 水量丰富, 易于造成突水灾害。由此分析, 本矿井充水的直接充水水源是主采煤层顶、底板砂岩裂隙水。在留设防水煤柱的情况下, 四含、五含水是矿井充水的补给水源, 而煤系下伏的石灰岩岩溶裂隙水则是矿井充水的重要隐患。

(二) 矿井充水通道

矿井充水的通道主要有:砂岩构造裂隙、采动冒落带裂隙、断层及岩溶陷落柱以及未封闭好的钻孔。一般来说, 地压水进入矿井的途径有以下几种:1) 高压灰岩水突破煤层底板进入矿井, 特别是在底板较薄的区段, 由于采动时底板造成破坏, 底板的有效隔水层厚度可能不足以抵抗下伏的高压水, 因此可能导致突水;2) 地下水 (包括老塘水) 沿断裂带进入矿井;3) 采、掘揭露或接近含、导水岩溶陷落柱导致突水;4) 沿封闭不良的钻孔沟通其它含水层导致突水。

四、结论

通过对含水层的充水因素分析, 得出不同含水层对矿井的影响程度不同。主要表现在:新生界松散层富水性中等, 只对浅部煤组的开采有影响, 容易控制;而各煤层顶、底板砂岩裂隙含水层一般处于封闭-半封闭状态, 出水量不大, 对矿井的开采影响的可能较小;煤系内断层在天然状态下, 处于相对稳定, 当随着开采的进行使各含水层之间发生水力联系时, 会导致突水的可能, 要加强观测和预防能力;矿区是逆断层发育, 断距较大, 使可采煤层与灰岩岩溶裂隙多处对口接触, 可能使矿井易发生较大突水事故, 应在开采前做好物化探工作, 保证突水事故在可控制范围内。

参考文献

[1]国家煤矿安全监察局.煤矿防治水规定释义[M].徐州:中国矿业大学出版社, 2009.

[2]杨坤光, 袁曼明.地质学基础[M].北京:中国地质大学出版社, 2009.

[3]施龙青, 韩进.开采煤层底板"四带"划分理论与实践[J].中国矿业大学学报, 2005.

[4]虎维岳.矿山水害防治理论与方法[M].北京:煤炭工业出版社, 2005.

含水层水文地质 第6篇

枯期径流是流域水资源重要组成部分,也是维持河流生态环境健康的重要因素。以饱和含水层地下水出流为主要成分的枯季径流特征受流域特征和含水层水文地质特性影响。利用流域出口断面枯季径流衰退过程可以推求流域水文地质参数。

Brutsaert和Nieber[1]通过对流域概化,根据枯季径流衰退特征和Boussinesq方程的解析解,将流域水文地质参数与退水过程建立了关系,成功地推求了6个流域的水文地质参数。Szilagyi等人[2](1998)运用数值模型验证了Brutsaert-Nieber方法确定水文地质参数的有效性。Brutsaert和Nieber提出的退水分析方法克服了实验室实验(如Darcy实验)与野外试验(如抽水试验)仅能研究厘米至几十米尺度上土壤及含水层特性的不足,将样本尺度的水文地质特性与流域、区域等较大尺度的特征建立了关系。该方法得到广泛应用,如Brutsaert和Lopez[3](1998)推求了Washita流域的22个子流域的渗透系数和给水度,并讨论了流域渗透系数与给水度的尺度效应。Zhang等人[4](2009)将Brutsaert-Nieber方法应用到我国东江的13个子流域,估算渗透系数和给水度。张艳芳[5]等结合Mendoza等人[6](2003)提出的方法,将BrutsaertNieber方法用于推求喀斯特区域的水文地质参数。

许多学者对Brutsaert-Nieber方法进行了改进和发展。针对Brutsaert-Nieber方法中流量下包线确定缺乏客观标准问题,Troch和De Troch[7](1993)提出了采用去除5%和10%数据点方法来确定下包线。针对流域坡度的影响,Brutsaert[8](1994)提出了考虑山坡坡度时的解析解,并分析了其对流域水文地质参数的影响。针对蒸散发对退水过程的影响,Szilagyi[9](2007)、Wang和Cai[10,11]等分别提出了根据退水特征估算流域蒸散发的计算方法。一些学者还将Brutsaert-Nieber方法用于分割基流[12]、评价人类活动对径流影响[10,11]。

上述方法均建立在饱和渗透系数垂向一致的均质含水层上。对于天然流域,由于土壤成土特征、含水层荷载变化、植被根系等因素的影响,饱和渗透系数常呈现随深度递减的特征[13]。饱和渗透系数垂向变化对地下水出流过程有着重要的影响。当考虑饱和渗透系数垂向幂函数变化时,Rupp和Selker[14](2005)给出了含水层快速退水和慢速退水过程的解析解(以下简称Rupp-Selker方程),并用数值方法验证了解析解的精确性。

本文利用Brutsaert-Nieber方法和Rupp-Selker方程,对东江星丰流域逐日退水过程进行分析,推求其水文地质参数,并重点讨论饱和渗透系数垂向幂函数变化时对水文地质参数推求结果的影响。

1 研究方法

1.1 退水方程

对于具有水平不透水底板的含水层(如图1),潜水水流运动可由Boussinesq方程描述:

式中,h=h(x,t)是地下水水位,t是时间,K是饱和渗透系数,f是给水度,x为横坐标,表示在含水层剖面的位置。

在含水层垂向为均质或等效均质时(K垂向不变),当初始含水层饱和,地下水出流过程未影响到分水岭处(图1中x=B处,B为含水层宽度)地下水水位(等同于式(1)中x→∞),Polubarinova-Kochina(1962)给出了含水层的出流公式:

式中,q为含水层单宽流量,D为含水层厚度。式(2)就是快速退水公式。

当分水岭x=B处潜水水位开始下降时(如图1),式(2)不再适用。Boussinsq(1904)用不完全Beta函数的逆函数描述初始潜水水位,给出了解析解:

这就是慢速退水方程。

当饱和渗透系数在垂向上呈幂函数变化时,即:

式中,z为相对于含水层底板的高度,KD为含水层表层饱和渗透系数。图2描述了不同n值时K的垂向变化。当含水层饱和带厚度为h时,饱和带平均饱和渗透系数为:

当h=D时,可求得含水层平均饱和渗透系数Kave:

Rupp和Selker[14](2005)给出了饱和渗透系数在垂向上呈幂函数变化时快速退水和慢速退水过程的解析解。对于快速退水:

式中:

式中:Beta(n+2,m+1)表示由n+2和m+1确定的Beta函数。m为无量纲参数,可使用K为常数时的q值校正m,Rupp和Selker[14](2005)建议m取值为1。

对于慢速退水:

其中BN是Beta函数:

当n=0,式(7)、(10)即为式(2)、(3)。

1.2 退水斜率曲线分析

天然河流径流消退过程中,流量Q是时间的函数,其变化率与流域调蓄特性有关:

当时间间隔为Δt,连续观测的流量数据为Qi和Qi+1,式(12)可近似写成:

式中,f()是反映流域调蓄特征的函数。相对于由降雨引起的其他径流成分,如河道直接径流、地表径流、大孔隙流等,地下径流的衰减速率|dQ/dt|是最小的。这就意味着在利用实测数据绘制图形时,可以采用取下包线的方法代表地下水退水过程。上述分析过程可称为退水斜率分析,其下包线称为退水斜率曲线。

退水斜率曲线f()呈现幂函数形式时,式(12)可写成:

式中,a和b是常数。

当流域概化成河道位于中间的矩形等效山坡时:

流域有效含水层宽度可表达为:

式中,A为流域面积,L为流域有效河长。

对于快速退水公式(7),对应式(14),a和b分别为:

式中:

a1与b1表示快速退水时的a、b值。对于慢速退水的式(10),a、b分别为:

式中:

a2与b2表示慢速退水时的a、b值。

1.3 水文地质参数反演

求得快速、慢速退水时的a、b值后,式(17)和(20)可联立方程组以推求水文地质参数。若已知流域基本特征,如面积A、有效河长L及平均含水层厚度D,即可反演流域饱和渗透系数KD和给水度f,计算公式为:

其中:

2 研究区域概况及退水过程分析

2.1 研究区域概况

星丰流域位于东江上游,流域面积42.6km2,流域平均宽度4.3km,平均比降3.14%,河流密度为13.62km/km2。流域气候属于湿润多雨的海洋性季风气候,多年平均降水量1720mm左右,年内降水主要集中于4~9月,降水年际、年内分配不均匀。星丰流域的土壤主要为麻赤红壤、页赤红壤和潴育水稻土,土壤类型概化为砂壤土、壤土和粘土3种类型,其中,砂壤土占71.17%,壤土占20.44%,粘土占8.39%。土壤厚度多在1.5~2.5m之间。

2.2 退水过程分析

本文使用星丰水文站1969~2008年40年逐日流量资料分析退水过程。退水过程选择应满足以下条件:(1)选择雨后尽可能短的时刻作为地下水退水的起始点[3],本文选择雨后第2日作为退水分析的起点;(2)在退水过程中,地下径流量和其变化速率(-dQ/dt)连续递减,应去除当日径流大于前日径流、当日退水速率大于前日退水速率的径流资料。经条件(1)、(2)处理后,在双对数坐标轴上绘制Q与-dQ/dt关系图,如图3(a)。可以看出在-dQ/dt较小时,数据点成水平层状分布,且非常分散,难以识别出慢速退水的下包线(图3(a)),下包线确定存在不唯一性。为此,Rupp和Selker[15](2006)通过把Δt作为一个变量实现了对流量数据低值部分的修正,计算方法如下:

同时,可使时间间隔ti-ti-j所对应的流量差值|Qi-Qi-j|大于某一阈值(本文阈值为0.05m3/s)。修正后Q与-dQ/dt关系图如图3(b)。

由于-dQ/dt有一定的误差,分散的数据点难以确定下包线的位置。Troch和De Troch[7](1993)提出以去除5%的数据点的方法确定下包线的位置,为检验这个位置所得结果的敏感性,同时考虑去除10%数据点确定的下包线位置。

经上述过程处理后,n=0(K垂向均质)、n=4(K垂向非均质性显著)时,去除5%和10%的数据点Q~-dQ/dt下包线如图4。

3 水文地质参数推求

本文首先分析n=0、n=4时,去除5%和10%数据点下包线位置以及含水层厚度变化对饱和渗透系数与给水度的影响。星丰流域土壤厚度在1.5~2.5m之间,取流域含水层平均厚度为1.5m、2m、2.5m时,所求KD、Kave与f值见表1。从表1可以看出:(1)随着D增加,KD、Kave和f均减小,当D从1.5m增加到2.5m,KD、Kave和f分别减小约64%、64%和40%;(2)表1中去除10%数据点所求KD、f值均小于去除5%数据点所求值;数据去除比例由5%变为10%时a1、n=0时的a2、n=4时的a2,分别增加45%、27%与24%,变幅较大,但KD、Kave值变化相对较小,约减小10%,f约减小为20%。

(线1、3、5和线2、4、6分别为去除5%和10%数据点所定的下包线)

注:表中括号内、外分别为去除10%、5%数据点确定下包线所反演得到的值。

如果流域含水层厚度取2m,下包线位置以去除10%数据点确定,以此分析考虑饱和渗透系数垂向变化对推求KD、Kave与f的影响。如表2,随着n的增加(饱和渗透系数垂向变化增大),KD、Kave与f均增大。n从0增大到4时,推求的表层渗透系数KD增大了15倍,含水层平均渗透系数Kave增大了2.2倍,给水度f增大了1.4倍,饱和渗透系数较给水度敏感的多。由表2和式(21)还可以看到,随着n的增加,含水层底部渗透性减弱,b2值增大,慢速退水曲线斜率增大。

n值变化情况下的水文地质参数变化

当D的变化范围在1.5~2.5m,以去除10%数据点确定下包线,KD与f随n的变化如图5。对KD值取对数后,对不同的n值,KD的变化幅度大体相等。f的变幅随着n的增加而增大,由n=0时的0.06增大到n=4时的0.15。随着n的增加,KD和f的上下限值也均增大。

星丰流域土壤类型主要为砂壤土(占71.17%),砂壤土的饱和渗透系数大约为1.411×10-5~4.234×10-5m/s[16]。当流域含水层厚度为2m,本文所推求的含水层平均饱和渗透系数(Kave)在此范围内(n<1)或略大于此范围(n≥1)。推求的给水度与Zhang等人[4](2009)在东江流域应用Brutsaert-Nieber方法所得结果相符,与Szilagyi[17](2004)、Hilberts等人[18](2005)的结果也较接近。

4 结论

本文通过对星丰流域地下水退水过程的分析,系统阐释了退水过程下包线确定方法及其对水文地质参数(饱和渗透系数和给水度)的影响。利用Rupp-Selker方程,揭示了饱和渗透系数垂向幂函数变化对Brutsaert-Nieber方法反演的水文地质参数的影响。研究表明:

(1)通过改变退水时段间隔Δt,可以减少确定Q~-dQ/dt下包线的不确定性。

(2)采用去除5%和10%数据点确定下包线位置来计算水文地质参数所产生的不确定性小于含水层厚度和垂向非均质性变化产生的不确定性。

(3)随着饱和渗透系数垂向非均质性的增加,慢速退水斜率增大,推求的给水度和饱和渗透系数均增大。

含水层水文地质 第7篇

贵刊2009年第11期刊载了肖有才教授等的“井损的消除方法及其在抽水试验中的作用”论文(以下简称“肖文”)[1]。肖文认为，在稳定井流中，对每一抽水降次分别用俩俩观测孔配对求算渗透系数k及影响半径R，并比较它们的相对误差，误差最小者即认为是试验观测精度最高者;当相对误差<2%时，这一抽水方法才称得上是“标准”方法，由此求算的井损值也最准确。其次，肖文认为，欲使抽水工作更节省财力，利用在消除井损后的主井与1个较远处的观测孔就可使求解参数的精度达到95%以上。

笔者认为，用俩俩观测孔配对解算的水文地质参数相对误差最小者来判定试验的观测精度有待商榷，用消除井损后的主井与1个观测孔资料配对求算的参数精度并不一定能达到95%以上。

在抽水试验中，水位、流量的测量误差一般难以避免。为尽可能地消除误差，利用全部孔的资料求算的水文地质参数才更具代表性、综合性(特别异常者除外)。当采用单孔求算参数时，根据观测孔抽水过程或恢复过程的s-t数据，使用最小二乘法求算水文地质参数应是更好的选择，它利用了更多的有用数据，可最大限度地消除水位测量误差。

为了便于讨论，有必要回顾一下裘布依井流模型及通用稳定井流模型。

1 稳定井流模型

1.1 裘布依井流模型

法国水力学家J·裘布依于1863年研究了一个均质各项同性的圆柱体含水层，其周边为圆形定水头的给水边界，抽水井(完整井)位于圆柱体中心，其含水层中的水位降深可由下式表达:

式中:Q———抽水井的稳定出水量(m3/d);

Rd———圆柱状含水层半径(m);

r———对观测井，为抽水井中心至观测井中心的水平距离(m);对抽水主井，为井半径(m);

M———含水层厚度(m);

K———含水层的渗透系数(m/d);

s———观测点的水位降深(m)。

这个模型中的Rd，即通常所称的“影响半径R”，实质上为圆形定水头供水半径。为区别起见，本文称之为“裘布依半径”。

1.2 通用无界稳定井流模型

自然界难以找到与裘布依井流模型一致的圆柱体含水层。为此，陈雨孙(1985)研究了具有补给的通用无界稳定井流模型，给出了在有补给(水平及垂向)的均质各向同性无界含水层完整井的抽水通用表达式[2]:

式中:I———地下水天然水力坡度;

θ———与地下水流向的水平夹角(°);

E———含水层的补给因子(m-1);

a———含水层的导水系数(m2/d);t———抽水延续时间(d)。

式(2)既可描述非稳定过程，也可描述稳定过程。当达到稳定时，可用下式表达:

式中K0(E·r)为0阶修正Bessel函数:

在给定允许相对误差ε≤2.18时，(3)式可近似表达为:

令

当天然地下水水力坡度较小时，(5)式就可写成:

对于抽水主井，式(7)中的r需用主井等效半径rwe代替[3,4]:

式中:rw为主井半径，rn为泥浆堵塞层的外半径，k为含水层渗透系数，kn为泥浆堵塞层的渗透系数。特别地，当kn=k时，rwe=rw。

式(7)在形式上与裘布依公式(1)相同，但必须特别注意，它并不等价于裘布依井流模型。因为式(7)的适用范围并不是整个Ry，而是在0.2Ry的范围内;在此范围之外，不满足给定的允许误差。

Ry刻画了含水层补给能力的大小，陈雨孙称之为“引用补给半径”，它并不是可以实际观测到的物理半径，不等价于“裘布依半径Rd”，也不是习惯所称s=0所对应的“影响半径R”(该R实际上是“漏斗半径”)。笔者认为，不应混淆“裘布依半径Rd”、“引用补给半径Ry”及“漏斗半径R”三者之间的区别。

2 通用无界稳定井流水文地质参数的解算方法

井流模型的水文地质参数识别，最常用的做法就是将抽水试验的观测降深sg与井流模型的理论降深s用最小二乘法原理，使目标函数:

最小，从而解算出水文地质参数。该方法一般称为最小二乘法，本文称之为曲线拟合法，以区别于使用最小二乘原理的回归方法。

式(9)中的理论降深si可根据情况采用不同的计算公式，如式(2)、(3)、(5)、(7)。

目标函数与水文地质参数的关系一般为非线性关系，解算水文地质参数通常较为复杂。目前的数学工具可使曲线拟合法求解水文地质参数变得十分简单。刘大海(2001)使用Mathcad数学工具的Given-MinErr求解块解算非稳定泰斯井流、潜水博尔顿井流、有越流的汉图仕井流的水文地质参数，取得了十分满意的效果[5]。

3 计算引例及计算结果

3.1 计算引例

肖文引例,实际上是1964年陕西省咸阳市橡胶厂水源地勘察的一个抽水试验[2,3,6]，在文献[2]中作为一个无界越流含水层作了介绍。该抽水试验，试验层为第2含水层，含水层厚20m，天然地下水水力坡度为0.0006(坡度较小，忽略不计)。为便于讨论，重载其抽水稳定时的流量Q及降深实测资料，见表1(与肖文所载略有不同)，其s-Q曲线图及回归方程见图1。

3.2 确定修正降深

一般由于过滤器、填砾、泥浆堵塞以及观测误差等原因，当s与Q为非线性关系时，在抽水主井中存在水跃值，在观测井中可能存在水位滞后效应及观测误差等。因此，不仅应对抽水主井，也应对观测井建立s-Q的回归方程，分离出满足达西渗流的修正降深s′，见表2。

3.3 水文地质参数解算

本文用两种方法解算水文地质参数:曲线拟合法及s-lnr直线回归法。前者选择式(3)、(9)，用Mathcad的Given-MinErr求解块求算;后者选择式(7)、(9)，用Excel求算。解算结果见表3、表4。

4 水文地质参数的误差分析

抽水试验中，水位、流量的观测误差通常难以避免，从而影响水文地质参数的解算精度。观测误差主要受抽水试验方法、抽水设备、流量及水位等测量工具的影响。

4.1 流量误差对修正降深的误差分析

采用多降次抽水试验确定修正降深s′时，流量的测量误差将对修正降深产生误差ds′:

流量的测量误差一般比降深的测量误差大得多，其误差通常可达到1%～5%。据此计算，引例第1降次的主孔修正降深误差或可达到0.09～0.46m，如此大的误差是不可忽视的。因此，由流量测量误差引起的降深误差必将对参数解算产生误差。

4.2 测量误差对水文地质参数的误差分析

由式(7)，可导出由于流量、降深、及孔距测量误差对水文地质参数解算的相对误差。

4.2.1 单孔测量误差对解算参数的误差分析

从式(13)、(14)可以看出:渗透系数k及引用补给半径Ry的相对误差，对于观测孔主要由流量、降深的相对误差组成;对于抽水主井，井径的相对误差也不可忽视。

主井r通常不考虑泥浆泥皮堵塞效应而取进水段滤水管半径或基岩裸孔(钻头)半径rw，但据式(8)，当堵塞效应较大，或基岩的实际孔径比钻头孔径大得较多时，用rw代替rwe时将产生较大的误差。因此，即使主孔消除了紊流及三维流影响，也可能由于主孔孔径误差导致引用补给半径Ry产生较大误差。

4.2.2 俩俩孔测量误差对解算参数的误差分析

式中:Δs=s2-s1为孔2与孔1的降深增量，Δlnr=lnr2-lnr1为孔2与孔1的对数孔距增量。

从式(15)、(16)可以看出:渗透系数k的相对误差主要由流量、两孔降深增量的相对误差组成;引用补给半径Ry的相对误差，与降深增量的相对误差成正比，与降深增量成反比。

4.3 参数解算的实际相对误差

使用不同的观测资料及不同的解算方法，其解算结果存在一定差异。解算的参数1与参数2的实际相对误差可用下式确定:

4.3.1 实测降深与修正降深资料解算结果的相对误差

当s与Q为非线性时，一般存在紊流、三维流、井孔的泥浆堵塞效应及观测误差等，此时不仅需要对主孔、也需要对观测孔降深进行修正以分离出满足达西渗流的修正降深。用修正降深资料解算的参数值为真值，用实测降深资料解算的参数值为近似值，解算参数的相对误差见表5。总体而言，参数k的相对误差远小于参数Ry的相对误差。利用全部观测孔(B1,B4,B5)资料解算的参数相对误差较小，参数k的相对误差<10%，而参数Ry的相对误差则达到20%～30%;利用带主孔资料解算的参数相对误差较大，参数k的相对误差10%～20%，而参数Ry的相对误差则达到20%～60%乃至更高。

4.3.2 直线回归法与曲线拟合法解算结果的相对误差

曲线拟合法采用了通用无界井流模型的解析公式，解算的参数值为真值。s-lnr直线回归法，采用的是通用无界井流模型的近似公式，仅适应于r≤0.2Ry的范围内，解算的参数值为近似值;在其范围外，不满足给定的精度要求。

直线回归法与曲线拟合法解算结果的相对误差见表6。参数k的相对误差较小，一般在10%以内;参数Ry的误差则较大，可达到10%～50%左右。在s-lnr直线回归法的适用范围内(不带B5孔)解算参数的相对误差显著小于适用范围外(带B5孔)解算参数的相对误差。采用修正降深资料解算的参数，其相对误差得到显著改善。

4.3.3 带主孔与不带主孔资料解算结果的相对误差

由于抽水主孔的观测误差及孔径误差一般远大于观测孔，一般而言，不带主孔资料解算的参数值更为真实，而带主孔资料解算的参数值代表性较差。带主孔与不带主孔资料解算的水文地质参数，其相对误差见表7。参数k的相对误差在10%～30%左右，Ry的相对误差在20%～70%左右。

5 若干问题的讨论

从以上分析及参数解算结果对比中可以得到以下几点认识:

(1)选择通用无界稳定井流模型及曲线拟合法，是解算水文地质参数的通用方法及精算方法，目前的数学工具使其解算并不困难。

(2)由于抽水主井一定范围内存在紊流、三维流、井孔的泥浆泥皮堵塞效应以及水位、流量的观测误差等，宜对所有孔进行降深修正以确定符合达西渗流的修正降深。

(3)渗透系数k及引用补给半径Ry的相对误差，对于观测孔主要由流量、降深的相对误差组成;对于抽水主井，井径的相对误差也不可忽视。

(4)Ry的解算误差远大于k的解算误差。

(5)参数解算的s-lnr直线回归法，适用范围为r≤0.2Ry。在此范围外，其参数解算误差不满足给定的精度要求。实际解算结果表明，观测孔B5位于s-lnr直线回归法给定允许相对误差的适用范围之外。在适用范围内，s-lnr直线回归法的解算结果与曲线拟合法的解算结果十分接近(表3);在适用范围外，两者的解算结果差异较大(表4)，Ry参数尤为如此。

(6)为了最大限度地消除水位、流量的测量误差，使解算的参数更具代表性，应尽可能多地利用观测井资料及其修正降深，使用最小二乘原理求算水文地质参数。利用俩俩观测孔配对解算的水文地质参数，其相对误差最小并不表明该参数的精度最高。

(7)抽水主孔，即使消除了紊流、三维流的影响，也可能由于抽水主井的修正降深误差及泥皮堵塞效应导致水文地质参数的解算误差过大，对Ry参数尤为如此。因此利用消除井损后的主孔与1个观测孔进行参数解算，其精度并不一定能达到95%以上。

参考文献

[1]肖有才等.井损的消除方法及其在抽水试验中的作用[J].工程勘察,2009.11.

[2]陈雨孙.颜明志.抽水试验原理与参数测定[M].北京,水利电力出版社,1985.8.

[3]陈雨孙.单井水力学[M].北京,中国建筑出版社,1977.12.

[4]杨金忠等.地下水动力学原理[M].北京,科学出版社,2009.12.

[5]刘大海.非稳定井流模型水文地质参数识别的Mathcad方法[M].深圳市地质学会2001学术年会论文集,中国地质大学出版社,2002.7.

含水层水文地质 第8篇

在地下水流数值模拟中,准确刻画含水介质非均质性的空间分布特征是得到可靠的数值模拟结果的前提。然而由于自然界中水文地质环境形成机制的复杂性,使得我们难以精确掌握描述含水层变化的水文地质参数。采用钻井采样或抽水试验等方法往往受到经费和技术等各种因素限制,只能获取少量的实测数据,而无法充分了解整个区域的参数分布。寻求一种尽可能利用有限的实测资料对未知位置的含水层参数进行合理估值的方法,是目前大区域地下水流模拟中的关键问题(宋刚等,2005)。从地质统计方法开始发展以来,该方法在水文地质参数空间变异性等方面得到了广泛的应用。Freeze(1975)经由众多现场调查所做的统计分析发现了渗透系数大多具有对数正态分布的特性。Delhomme(1978)进一步研究发现:经过对数转换的水文地质参数,在正态分布下,具有较原值更佳的空间变异结构,故后来的学者大多采用对数形态数据,作为分析试验变差函数的基本数据。Delhomme(1979)、Hoeksema和Kitanidis(1985)、Woodbury和Sudicky(1991)、Eggleston等(1996)、Eggleston和Rojstaczer(1998)等先后以试验变差函数及Kriging法推估水文地质参数空间结构并评估其不确定性。

由于含水层的水文地质参数还常常受到水文地质环境形成时的影响,具有各向异性的特征(施小清等,2005),这更增加了掌握水文地质参数分布情况的难度。本文为解决水文地质参数空间结构不易获取的难点,依据上海市地面沉降项目收集的大量实测资料为基础,采用地质统计方法,将渗透系数的自然对数值作为探讨分析的区域化变量,深入分析研究区的水文地质空间结构分布与变异性,结合研究区的水文地质资料加以佐证,作为地下水流数值模拟的重要前提步骤。

1 研究区域概况

上海市位于长江三角洲东缘,属长江三角洲前缘河口滨海平原和太湖流域东缘低地,陆域总面积6430km2,水域面积约1700km2,其中市区面积370km2。区内除西南部残留有高出地面数十米至百米的零星孤丘外,一般地势平坦,略成东高西低的倾斜状平原。地面高程(吴凇口基准面)一般为2.2~4.8m。

钻孔资料表明(宋保平等,2004),研究区埋深150m以下地层为下更新统,岩性主要为杂色粘土和砂互层,属河流-湖泊沉积;埋深150m以上属于中更新世、晚更新世和全新世沉积,岩性主要为灰色粘土和灰色砂、粉砂相间,属海陆交互相沉积。上海地下水赋存条件和分布规律受控于区域地貌、地层岩性及厚度,根据研究区地下水系统特点,可将其划分为1个潜水含水层,5个承压含水层,6个弱透水层(隔水层)。

受资料限制,本次选取了研究区内通过抽水试验得到具有较多实测渗透系数观测值的第三承压含水层作为研究对象。该含水层实测数据共184个,观测点位置分布和数值如图1所示。由图中可以大致发现观测数据分布不均匀,主要集中于市区,从数据大小来看,总体上西南区渗透系数值较小,而东北部较大。

2 观测数据分析

由于研究区内部分实测渗透系数值较低,且最大值与最小值间相差百倍至千倍(见图2a),为防止在取对数后出现负值,首先将渗透系数单位由m/min统一转换为m/d。Delhomme(1978)建议通过自然对数转换含水层的渗透系数后,可得到较佳的空间变异结构,因此对渗透系数K取对数后(ln K)进行分析。表1为该含水层渗透系数观测值基本统计结果。统计结果显示,渗透系数的范围在183.19m/d至1.63 m/d之间,偏态为正偏;取对数后范围在5.211 ln(m/d)至0.489 ln(m/d)之间,平均值约为2.781 ln(m/d),差异不大,且偏态为负偏。将渗透系数取对数变换前后的频率直方图进行对比(如图2所示),同时将转换前、后的渗透系数进行卡-方检定(Chi-square test),并由组数及参数的个数计算条件下,检验正态分布优度,经计算后可知取自然对数后的渗透系数自由度(d.f.)为5,其统计值为3.531(见表2),在5%显著水平下其理论卡方检定值为11.0705,由此可知该含水层的渗透系数呈对数正态分布,不满足正态分布。

3 空间变异结构分析

3.1 试验变差函数(Experimental Variogram)的计算和分析

本研究选取GSLIB(Deutsch和Journel,1998)作为分析工具,首先获取全方位试验变差函数(omni-directional variogram)作初步变异结构和理论模型的研究,即包含所有观测数据进行计算,以获取一个较为平均的结果,也可以从图形中判定渗透系数是否符合平稳假设。

由于观测数据之间的距离一般均大于2.5km,选取间距(Lag distance)为2.5km的情况下计算该含水层对数渗透系数的试验变差函数,并将间距(h)与计算出的试验变差函数γ(h)绘成全方位试验半变异图(图3)。由该图可以发现观测点在短距离内即有很大的变异产生,且在末端有部分数据点呈现上下跳动的现象,理论上变差函数接近基台值时,应有趋于平缓的趋势,造成此现象的原因,可能是选取的样本数太少,以至于数据对较少的间距范围内易受局部变化较剧烈的数值影响,即受到极值的影响,而产生较高的不确定性。

一般水文地质上常采用球状模型和指数模型进行试验变差函数的拟合(王仁铎和胡光道,1988;朝伦巴根等,1994),本次同时采用上述两种模型对全方位试验变差函数进行拟合,结果如图3和表3。结果表明:两者拟合结果相差不大,研究区变程约为17~18km。考虑到本例指数函数的决定系数(coefficient of determination)略高,故下文均采用指数函数进行试验变差函数拟合。

注:基台值等于块金效应和拱高之和。

3.2 各向异性分析

当区域化变量在空间的变异性仅随着距离而改变,而与方向无关时,称为各向同性,否则则称为各向异性。由于受观测点分布位置的影响,在针对某一特定方向上计算试验变差函数可能出现数据对数不足的情况,因此实际中的做法是在采用GSLIB软件的计算时,选定一容忍角(Tolerance Angle),而后于特定方向上加减此容忍角,并在此范围内按照变差函数公式加以计算,以获得该选定方向的试验变差函数,并用理论变差函数模型拟合(Deutsch和Journel,1998)。

在分析各向异性分析中,需要确定两个主要参数:

(1)容忍角的确定

若容忍角越小,则该选定方向上计算所得到的试验变差函数越精确。但容忍角过小则包含的数据配对数太少,分析结果往往表现为不稳定。最佳的方法是尝试不同大小的容忍角,选择其中使得试验变差函数结果较为稳定的最小容忍角。

(2)距离参数的确定

在计算试验变差函数时候,需要确定两个距离参数,即滞后间距(lag increment)和容忍间距(lag tolerance)。若取样点分布均匀,则取样点间的间距即为最合适的分级间距,否则,则应以观测取样点间的平均距离作为初始的分级间距,再根据取样点的多少,调整间距大小至得到合理、稳定的试验变差函数图。为获得较稳定的试验变差函数图,各间距所使用的数据对数应超过30对以上(Journel和Huijbregts,1978;王仁铎和胡光道,1988)。本研究所使用间距根据取样数取为2500m,使得各间距所使用的对数均超过30对。容忍间距则取为分级间距的一半。

在选定容忍角、距离参数后,需要分析确定各向异性主轴,主要步骤如下:

(1)计算多个不同方向的试验变差函数;

(2)拟合试验变差函数,获得理论变差函数模型在各个方向的变程;

(3)根据上一步骤结果,绘制玫瑰图(Rose Diagram)。若区域化变量具有各向异性时,玫瑰图呈现椭圆形,其长轴方向表示最大连续方向,短轴表示最小连续方向。当有套合结构(即不同地质类型)时,玫瑰图的长短轴大小及方向均不同。

为了检验研究区对数型渗透系数空间变异结构是否具有各向异性,本次计算了多个不同方向的试验变差函数;选取了正东向、北偏东60°、北偏东30°、正北向、北偏西30°、北偏西45°、北偏西60°、北偏西75°8个方向作为分析方向。在计算过程中选择尝试了不同容忍角11.25°和22.5°检验。但由于11.25°包含的观测数据对太少,以至于无法获得稳定的结果,因此容忍角改为22.5°,计算结果见表4。绘制由各个方向的拟合模型所得变程的玫瑰图(图4),由图中可以看出,最小变异方向为北偏东30°,最大变异方向为与其正交的东偏南30°。

分析结果同时显示,在北偏东0~90°范围内的各个方向具有相近的基台值,但却有不同的变程范围,这是几何各向异性的特征;同时东偏南0~90°范围内的各个方向也与此相似。

从地质角度来看,研究资料表明(陈颐亨,1988;宋保平等,2004):上海地区正处于苏北凹陷的南缘,所以本区的第四系基底地形表现为西南高、东北低。这一基底地形奠定了本区古河道和三角洲展布和沉积物发育的基本格局,也确立了本区地下水动力场的宏观格局。同时东偏南30°为本区地层的沉积方向和冲淤积方向,因而造成该方向变异程度大,而与其垂直的南北方向变异程度小。故而渗透系数整体上表现为由西南部湖沼低地富水性较弱向东北方向渐趋增强的趋势,这正说明了研究区各向异性方向与地质环境形成过程息息相关。

3.3 交互验证

在进行克里金插值推估前,需先针对前文选定的理论变差函数模型进行验证,检验原始变差函数模型套用克里金交互验证(Cross-Validation)所产生的残差值(Reduced Residuals)(原始数据减去插值)的概率密度函数是否满足标准正态分布,即克里金平均误差(Kriging Average Error;KAE)是否接近零,克里金均方误差(Kriging Reduced Mean Square Error;KRMSE)是否接近1作为验证标准(张仁铎,2005),若方差大于1,表示实际方差比估计方差大,需将变差函数的基台值调高做校正;反之,需将基台值调低,重新验证残差值的统计分布是否接近标准正态。交互验证程序使用GSLIB中的KT3D程序进行运算,过程中输入前文变差函数分析所得的空间信息的长短轴方位与长度等数据进行交互验证计算。交互验证的最终结果KAE值为0.032,KRMSE值为1.013,皆在允许误差范围之内(Chiles和Delfiner,1999)。

4 Kriging插值

经交叉验证检定理论变差模型的可行性后,在研究区域内采用KT3D程序以95行97列(0.9 km×0.9 km)的格网对第三承压含水层渗透系数的对数值进行二维普通克里金插值,结果如图5所示。可以看出,在研究区东部、东北部数值较高,中部大部分较低,西南部最低。这与研究区第三承压含水层的实际水文地质条件也较为符合。为了进一步检验克里金插值的效果,将插值后所得到所有网格点的数值和实测观测值进行统计分析对比,描述统计量对比见表5。

据表5显示可以看出,普通克里金插值结果较原始观测数据分布范围更小,标准差偏小,存在着平滑效应,即低值高估,高值低估现象。这也可以从频率分布的累计概率分布函数(图6)可以看出,虽然普通克里金插值结果的统计特性在整体变化上与实测数据分布较为接近,但插值结果表现为大部分集中于中间值,取高值和低值的概率很低,这正说明了以线性最优无偏估计为准则的普通克里金插值对于参数的空间分布插值较为平滑,更适用于空间变异性较小的区域,对于极值的推估能力较差。

5 结论和建议

本文的空间结构分析结果表明,上海地区第三承压含水层的渗透系数经过自然对数变换后,表现为正态分布,在平稳区域化变量假设的条件下具有显著的各向异性。在东偏南30°方向具有较剧烈的变化,而在与其正交的北偏东30°方向则变化较为平缓。研究区各向异性的方向与地质环境形成过程息息相关,西南高、东北低的第四系基底地形决定了本区渗透系数由西南部较弱而向东北方向渐趋增强之势。

对选定的理论变差函数模型进行交互验证后,采用普通Kriging法推估研究区的渗透系数分布,合理的反应了实际水文地质条件,弥补了观测资料缺乏的不足。同时研究结果发现普通克里金插值结果存在着平滑效应,即低值高估,高值低估现象,说明了以线性最优无偏估计为准则的普通克里金插值对于参数的空间分布插值较为平滑,更适用于空间变异性较小的区域,对于极值的推估能力较差。

需要注意的是,在渗透系数空间结构分析时,为克服某些方向上缺乏足够的观测数据造成数据对稀少的缺点,采用了较大的容忍角度进行分析,虽获得较稳定的变异结构,但容易受到其它方向的配对而混淆该方向的变异性,同时造成试验变差函数出现振荡,所以,建议将来可在数据对稀少的方向结合克里金估计方差寻找最佳的抽水试验井设置地点,以求获得更多位置的渗透系数值,从而在异向性分析时能获得更精确的结果。

摘要：本文采用地质统计方法分析上海第三承压含水层渗透系数的空间变异性。将取对数后的渗透系数作为区域化变量,计算各个方向的试验变差函数,进行理论变差函数模型拟合,以分析含水层空间结构并判定含水层可能的各向异性。变差函数模型经交互验证后,再采用普通Kriging法进行空间插值获得该含水层的渗透系数空间分布。研究结果表明:对数渗透系数在空间上具有明显的几何各向异性,大约在北偏东30°方向有最小的变异性,而在东偏南30°方向具有最大的变异性,这与含水层形成的水文地质环境具有密切关系。