Kriging插值模型在地下水位监测网优化中的应用

第 41 卷 第 9 期 2 0 1 0年 5月
人 民 长 江 Yangtze R iver
Vol . 41, No. 9 M ay, 2010
文章编号 : 1001 - 4179 ( 2010 ) 09 - 0014 - 04
Kriging插值模型在地下水位监测网优化中的应用
3
1. 5 临界理论方差 σ0 的确定
分析地下水监测网密度时 ,临界方差值的选取 ,目 前我国尚无相应规范规定 , 就实际情况以及国内外文 2 献来看 ,一般当理论方差 σ0 取 0. 5 ～0. 6 时 ,地下水监 测井网精度就可满足实际需要 。
2
1. 3 阶数 N 的确定
理论上计算某一估计值时 , 整个研究区的水位都 与该点水位相关 ,但实际上只有靠近估计点“ 内圈层 ” 的若干监测井才与估计值密切相关 , 而超过一定距离 的“ 外圈层 ” 井不但没什么影响 ,而且容易产生病态方 程组 ,而使计算失真 。为解决这种“ 屏蔽 ” 作用 , 实际 研究中宜采用最近的内圈层 2 ～6 个监测井参与计算 即 N = 2 ～6,并适时调整参与计算 N 值 , 使效果达到 最优 。
2 黄水河流域地下水监测网优化 2. 1 流域概况
黄水河流域位于山东半岛北部 ,地形是东南高 、 西 [3] 北低 , 南 部 为 低 山 丘 陵 , 北 部 为 平 原 , 流 域 面 积 2 1 034. 47 km ,有地下水水位长期监测井 22 眼 。 区内主要含水层为第四系孔隙潜水含水层 , 主要 分布于黄水河冲洪积平原区 , 以大气降水及河川径流 入渗补给为主 ,东部 、 南部山丘区基岩地下水侧向补给 量较少 。 由于受大气降水季节 、 年际分配极为不均的影响 , 本区地下水动态变化较大 , 年内枯水季节地下水位下 降 , 丰水季节地下水位上升 ; 年际连续枯水年地下水 位大幅持续下降 , 连续丰水年 , 地下水位回升幅度较 大。
© 1994-2010 China Academic Journal Electronic Publishing House. All rights reserved.
16
人 民 长 江
2010 年
( 1 ) 变差函数的确定 。首先计算两两井距 hij 并
监测地下水位理论方差整体偏大 ,而局部偏小 ,监测网 密度整体较低 ,而局部偏大 ,尤其是以地下水库边界为 界限更为明显 ,界限以外监测网密度偏小 ,界限以内监 测网密度偏大 。 根据区域优化目标 ,结合该区域的地下水动态 ,设 计备选方案为 : ( 1 ) 基本方案 。在密度偏低的位置增加监测井
刘 治 政 ,吴 晓 东 ,林 洪 孝
(山东农业大学 水利土木工程学院 ,山东 泰安 271018 )
摘要 : 目前运行的地下水监测网存在着很多问题 ,很难适应现代社会对信息的需求 ,需要对地下水监测网进行 优化 。详细介绍了 Kriging插值模型的建立及求解过程 ,并运用该模型对黄水河流域平原区进行了地下水水 位监测网的优化 ,得到了备选方案 ,然后根据具体的优化目标选出最优方案 ,结果较为合理 。 Kriging插值模 型充分考虑了水文地质条件的影响 ,是优化地下水监测网密度的有效方法 。 关 键 词 : 地下水监测网 ; 监测井 ; Kriging插值模型 ; 理论方差 ; 黃水河流域 ; 中图法分类号 : P641 文献标志码 : A
收稿日期 : 2009 - 12 - 24 基金项目 : 科技部重大国际合作项目 (2007DFB70200)
式 ( 2 ) 用来求 λi , 在 γ( xi , xj ) 已知条件下 , 是一 正定方程组 ,因此有唯一解 。 利用式 ( 2 )并注意方差函数 、 协方差函数 、 变差函 数在本征条件下的转换关系 , 可得计算误差的理论方
Z
3
[2]
( x0 ) =
6 λZ (x )
i i i =1
N
( 1)
3 式中 , Z ( x0 ) 为 x0 的估计值 ; λ为 Kriging权系数 。 利用式 ( 1 ) ,在保证无偏性和最优性的前提下 , 结
合协方差定义 ,并引入拉格朗日算法 ,可得 :
6 λγ( x , x )
i i j i =1 N
N
+μ = γ( xi , x0 ) ( 2)
6λ
i =1
i
=1
式中 , γ( xi , xj ) =
1 V a r [ Z ( xi - xj ) ] 为变差函数 ; Var 2
[ ]为方差算符 ; μ为拉格朗日算子 。
1 Kriging插值模型 1. 1 原理
Kriging方法是一种对时空分布变量求最优 、 线
测值可求出实验变差函数 γ ( h ) 为 : N ( h) 1 3 2 γ ( h) = [ Z ( xi ) - Z ( xi + h ) ] ( 4 ) 2 N ( h) 6 i =1 3 利用式 ( 4 ) 求出不同间距 h 所对应的 γ ( h ) , 再作 3 γ ( h ) 关于 h 的曲线拟合 , 从而得到 γ( h ) 。理论上 γ( h ) 曲线服从幂函数 、 高次多项式 、 高斯函数 、 球状函 数等分布形式 。 变差函数 γ( h ) 是刻画 Z ( x ) 空间上变化的规律 性的单调递增函数 ,它与地下水流系统所处的水文地 质条件密切相关 , 在计算 γ( h ) 过程中要注意与有关 条件相结合 。
0. 013 1 hij
1. 1826 2
,相应的变差函数曲线见图 1。
填 14 号井 ; ( 4 ) 方案 3。在基本方案基础上 , 封填 9 号井和 14 号井 。 为选出最优方案 , 需计算基本方案下各井的理论 方差 ,方法同上 ,部分计算成果见表 2。 从计算结果中可以看出 ,在 14 号井没有参与计算 的情况下 ,各相关井的理论方差仍然可以满足精度要 求 ,根据优化目标要求 , 14 号井可以封填 。 是否封填 9 号井 , 主要取决于封填前后与 9 号井 有关的监测井的理论方差是否满足精度要求 , 封填前 后部分对比结果见表 3。 如表 3 所示 ,保留与封填 9 号井对相关井的理论 方差是有影响的 ,虽然保留 9 号井会使 11 号井理论方 差偏低 ,但封填后井不能满足精度要求 , 根据优化目 标 , 9 号井应该保留 。综合分析 , 保留 9 号和 14 号井 虽然能够使监测网达到更高精度的检测标准 , 但监测 费用也相应增加 ,而将 9 号和 14 号井全部封填后 , 虽 然检测费用降低了 ,但监测精度又达不到优化目标的 要求。 综合考虑当地实际情况 , 宜采用方案 2 , 即在基
地下水监测资料是揭示地下水运动规律的重要基 础 ,目前全国各地已经建成相当规模的监测网 ,但随着 水资源开发利用和环境保护矛盾日益突出 , 监测网所 具有的弱点逐渐暴露 。 ( 1 )监测井空间布局的不合理性 。监测井空间位 置人为的随机确定 ,使得监测网没有明确的系统框架 , 监测井的空间布局缺乏合理选择 , 导致许多监测井的 空间位置缺乏稳定性 、 时间上缺乏延续性等问题 。 ( 2 )监测井数量的不合理性 。最典型的表现为部 分地段监测网井过密而造成数据信息的大量冗余 , 而 偏离城区特别是交通不发达的地区 , 监测井往往又控 制性不够 。 ( 3 )监测井网不能应对地下水动态的影响 。主要 表现为地下水位的持续下降 , 使得部分监测井失去监 测功能而被淘汰 。 因此 ,依靠目前的地下水监测网不可能获得可靠 、 符合要求的数据信息 ,需要对现行的监测网进行优化 , 以获得包含更多有效信息的监测数据 。
作者简介 : 刘治政 ,男 ,硕士研究生 ,主要从事水资源开发利用与保护管理研究 。 E - mail: liuzhizheng1985@126. com 通讯作者 : 林洪孝 ,男 ,教授 ,博士 ,主要从事水资源开发利用与保护管理研究 。 © 1994-2010 China Academic Journal Electronic Publishing House. All rights reserved.
1. 2 变差函数
变差函数是刻画变量 Z ( x ) 空间统计结构的 , 它 只依赖于空间点间的相对位置与随机场特性 , 这是 Kriging方法的理论基础 。 在 Z ( x ) 服从本征条件时 ,理论变差函数 γ( h ) =
1 2 E [ Z ( x ) - Z ( x + h ) ] 。当两点间相对距离为 h 的 2 实测数据有 N ( h ) 对时 ,利用 xi 点与 xi + h点上变量实
第 7期
刘治政 ,等 : Kriging插值模型在地下水位监测网优化中的应用
15
差σ 为: σ2 =
2
γ( x1 , x0 )
λ1 λ = , λ2 ┇ λ N
, E =
1 1
6 λγ( x , x )
i i
N
0
+μ
2 0
( 3)
i =1
B =
γ( x2 , x0 ) ┇ γ( xN , x0 )
1. 4 Kriging方程组的求解
为方便求解式 ( 2 ) ,在变差函数已知的情况下 , 将 其转化为矩阵形式 ,结果如下 : ( 5) Aλ +μE = B 式中 , γ( x1 , x1 )
A =
γ( x1 , x2 ) γ( x2 , x2 ) ┇
… γ( x1 , xN ) … γ( x2 , xN ) … ┇ … γ( xN , xN )
23 眼 。 ( 2 ) 方案 1。在基本方案基础上 , 保留 9 号井和 14 号井 ; ( 3 ) 方案 2。在基本方案基础上 , 保留 9 号井 , 封
按距离间隔将其分 5 类 ,即 hij ≤30 km , 30 km < hij ≤ 60 km , 60 km < hij ≤90 km , 90 km < hij ≤120 km , hij > 120 km。然后计算出每一类中 hij 的平均值以及区 域化变量 zi 、 zj 的平方增量 ( zi - zj ) 的平均值 。把 hij 均 1 2 ( zi - zj ) 均 作为 h 和 γ( h ) 的对应值 ,点绘实际曲 和 2 线 ,根据这些点在曲线上的分布 ,选取合适的理论变差 函数模型 ,用最佳拟合技术确定拟合曲线 ,从而得到合 理的变差函数 γ( h ) 。最终确定变差函数为 γ( hij ) =