地质统计学理论在水文数值模拟中的应用

发布网友发布时间：2022-04-28 11:05

共1个回答

热心网友时间：2023-10-02 05:15

邢永强¹　李金荣²　杨振放²

（1.河南省国土资源科学研究院，郑州 450016；2.郑州大学环境与水利学院，郑州 450001）

《安徽农业科学》，文章编号：0517-6611-(2007)-14-04101-02

摘要在地下水资源评价中，用数值方法进行水流模拟时，需要给出每个节点上含水层底板标高值。将地质统计学应用于空间分析学科，已被越来越多的学者所利用。文中着重阐述了地质统计学的基本原理及它们在含水层底板标高估值中的应用，通过计算结果可以得出克立格方法是进行含水层底板标高估值的空间最优估计方法。

关键词数值模拟区域化变量变差函数克立格方法

人类的各种活动使地表水或地下水的水质和水量发生很大的变化，因此为保护水环境我们需要进行水资源评价。水流模拟是水资源评价的一个核心内容，往往采用有限元或有限差分法（李俊亭，1989）进行评价，但是由于实际条件的*，我们掌握的节点数据往往是有限的，因此需要根据已知点的数据来求未知点的资料信息。传统的插值方法是根据已知测量值按人工线性插值查出其他节点上的变量值，其速度慢且精度低，况且没有考虑研究对象的空间变异性，而地质统计学中的克立格方法正好克服了这一缺点（孙洪泉，1990）。

克立格方法最初是由南非矿山地质工程师D.G.Krige（克立格）于1951年提出并以他的名字命名的方法，随后由法国学者完善并发展形成目前的地质统计学理论（孙洪泉，1990）。地质统计学最早应用于脉状矿山品位、储量研究；随着研究的深入，该方法广泛应用于矿产资源评价、钻探工作、采样工作；另外，还应用于化探、冶金、土壤（李恩羊等，1989）等研究领域。20世纪80年代之后，克立格方法在水文地质学领域得到了广泛应用（许多项，1993）。李金荣等（2003，2002）采用克立格法对含水层底板标高进行估值，得出克立格方法是对空间变量进行估值的一种较好的方法。

1 地质统计学的基本理论

1.1 区域化变量

区域化变量理论是地质统计学的理论基础。它的定义（孙洪泉，1990）是：以空间点x的3个直角坐标x_u，x_v，x_w为自变量的随机场Z（x_u，x_v，x_w）=Z（x）称为一个区域化变量。区域化变量在观测前，可以看做是随机场；观测后就得到Z（x）的一个实现，每一个实现Z（x）就是一个普通的三元实函数（或空间点函数）。它的显著特征是具有随机性、结构性。

含水层底板标高可以看成是区域化变量，它具有随机性和结构性。所谓随机性是指空间点x固定后，Z（x）为一随机变量，所谓结构性是指不同点x与x+h 处的 Z（x）与Z（x+h）之间具有空间相关性。

1.2 变差函数

变差函数是克立格方法计算的基础。它描述区域化变量的空间结构性，也描述其随机性。我们把区域化变量Z（x）在x和x+h两点处样本值之差的方差之半定义为Z（x）在x方向上的变差函数（又称变异函数），记为γ（x，h），即

环境·生态·水文·岩土：理论探讨与应用实践

实际工作中我们用实验变差函数γ^*（h）的公式计算：

环境·生态·水文·岩土：理论探讨与应用实践

式（2）中N（h）为相对于空间向量h、区域化变量Z（x）变异的统计点对数。对于不同的h可计算出相应的γ^*（h）值。

为了对区域化变量的未知值作出估计，还需要将实验变差函数拟合成相应的理论变差函数模型。拟合的模型有多种，其中常用的模型为球状模型（孙洪泉，1990）：

环境·生态·水文·岩土：理论探讨与应用实践

式中：C₀为块金常数，由微观结构变化和观测误差所决定的一种随机变化成分；C为拱高，Z（x）结构变化的最大值；C₀+C为基台值，反映一定方向上Z（x）结构变化与随机变化的总的变化幅度；a为变程，反映Z（x）变化的影响范围或变异速度。

1.3 克立格方法

从数学上讲，克立格法是一种对空间变量分布数据求最优、线形、无偏内插估计量的方法。它是以反映变量空间结构特征的变异函数为基础，以取得估计方差最小为目标，在无偏性约束条件下求优的一种估计方法。无偏性消除了系统误差，估计方差最小则表明有较高的精度，所以称最优估计。从水文地质角度讲，它根据已知观测点上水文地质变量（本文中针对含水层底板标高）的实测数据，对水文地质变量进行结构性（变差函数的模型确定）分析之后，为了对待估点作出一种线形、无偏、最小方差的估计，而对周围已知点的测量值赋予一定的权系数，进行加权平均来估计待估点水文地质变量的方法。

令区域化变量底板标高为B，设位于x₀点处含水层底板标高变量的待估值B^*（x₀）是周围全部已知的测量值B（x_i）（i=1，2，…，n）的线形组合，则：

环境·生态·水文·岩土：理论探讨与应用实践

式中：λ_i为克立格权系数（i=1，2，…，n）；n 为已知观测点总数。

式（4）中待定权系数λ_i（i=1，2，…，n）的确定应满足无偏性条件和最优性条件，才能保证估计量B^*（x_i）的线形、无偏、最优估计。

假定区域化变量底板标高B满足本征假设（孙洪泉，1990），即

环境·生态·水文·岩土：理论探讨与应用实践

（1）根据无偏性条件

环境·生态·水文·岩土：理论探讨与应用实践

经推导可得

环境·生态·水文·岩土：理论探讨与应用实践

（2）根据最优性条件（即估计方差最小条件）

环境·生态·水文·岩土：理论探讨与应用实践

欲使得式（8）估计方差

在满足无偏性条件下达到极小，就要利用求条件极值的拉格朗日乘数法，即要求：

环境·生态·水文·岩土：理论探讨与应用实践

达到极小时的诸λ_i（i=1，2，…，n）和μ。式（9）中F是n个权系数λ_i（i=1，2，…，n）和μ的（n+1）个元函数；μ 是拉格朗日乘数。

令

环境·生态·水文·岩土：理论探讨与应用实践

由此推导可得普通克立格方程组：

环境·生态·水文·岩土：理论探讨与应用实践

式（11）共包含了n+1个方程，可求解出n+1个未知的λ_i（i=1，2，…，n）和μ。

根据求出的λ_i（i=1，2，…，n）和μ，则可得到普通克立格方差表达式：

环境·生态·水文·岩土：理论探讨与应用实践

值得注意的是，式（11）和式（12）均只取决于变差函数及信息点与待估点之间的相对几何特征，而与信息样品数据B（x_i）（i=1，2，…，n）无直接关系，故它能够预测估计精度。

最后根据求出的λ_i（i=1，2，…，n）和μ，根据式（4）可计算出任一点的含水层底板标高值B，根据式（12）可得出任一点的估计方差。

2 应用

以陕西省咸阳某一水源地地下水资源评价为例。该研究区处于冲积平原上，评价的目的层为第四系含水层。用有限元法对该区地下水进行数值模拟，采用不规则三角形网格剖分，将全区剖分成674个单元，共382个节点，总面积为203km²（图1）。根据野外实测资料，取得含水层的底板标高系列数据n=68，作为含水层底板标高等值线图的绘制依据，利用前面所述的地质统计学方法理论，编制了相应的计算机程序，用该方法理论可计算出各节点的含水层底板标高估值和克立格估计方差，绘制了含水层底板标高等值线图（图2）。

图1 研究区单元剖分图

Fig.1 The unit analysis map of research area

从计算模拟结果（由于计算较多，这里不再一一列出）可以看出，用地质统计学方法计算出来的底板标高值与其实测值误差较小，最大达10.82m，最小为0.07m，平均为3.47m。用克立格方法计算的各实测点上的估计方差较小，最大达1.394 2，最小为0.003，平均为0.301 6。其绝对误差和相对误差值也很小。因此说克立格方法是一种较好的空间模拟方法。

图2 用克立格方法模拟的底板标高等值线图（单位：m）

Fig.2 The isoline map of lower bed level estimation simulated by Kriging method

3 结论

（1）区域化变量理论能够反映变量空间分布的随机性与结构性特征，变差函数是地质统计学的基本工具，是对变量空间分布特征估计最为有效的方法之一。从模拟计算的结果可以得出克立格法是一种最优、线性、无偏的模拟空间变量的方法。

（2）含水层底板标高这一变量具有明显的空间变异性。对这种变量的最优估计方法是以能反映含水层底板标高空间分布的随机性与结构性特征的克立格方法的基本理论为基础，在无偏性约束条件下寻求估计方差最小的一种克立格方法。

（3）由于克立格方法有不依赖于信息样品数据（含水层底板标高）的特征，在事先不知道待估点实际值的情况下，直接给出估计精度，且估计精度较高。

参考文献

李俊亭.1989.地下水数值模拟.北京：地质出版社.

孙洪泉.1990.地质统计学及其应用.北京：中国矿业大学出版社.

李恩羊，袁新.1989.作物需水量的最优估计.水利学报，（10）：45～49.

许多项.1993.数学插值方法在水文地质学中的应用探讨.长春地质学院学报，（10）：423～429.

李金荣，杨振放，李云峰，李金玲.2003.两种方法在地下水位估值中的应用.水文地质工程地质，（3）：42～46.

李金荣，杨振放，郭建青.2002.变差函数在地下水位估值中的应用研究.西北水资源与水工程，13（4）：6～9.

Utilization of Geostatistics in Hydrologic Numerical Simulation

Xing Yong-qiang¹　Li Jin-rong²　Yang Zhen-fang²

（1.Sciencial Researchinstitute of land and resource of Henan Province，Zhengzhou 450016；2.College of Environment and Water Conservancy，Zhengzhou Univ.，Zhengzhou 450001）

Abstract：In many cases of water resource assessment.，When numerical method is oftenadopted in groundwater flow modeling，lower bed level value of aquifer on each node should be given.The geostatistics is used in spatial analysis，it has been advocated by more and more experts.The geostatistics method and their basic application in estimating lower bed level value of aquifer is expounded.Through the result，the author can draw an important conclusion：Kriging method is spatial optimal estimation one in estimating lower bed level of aquifer.

Key words：Numerical Simulation；regionalized variable；variogram；Kriging method