各向异性.doc
五、各向异性Anisotropy 当区域化变量在不同方向呈现不同特征时,半变异函数在不同方向也具有不同的特性。我们称这种现象为各向异性。常见的各向异性有两种。 1几何各向异性Geometric Anisotropy几何各向异性的特点是半变异函数的槛值不变,变程随方向变化。如果求出任一平面内所有方向上的半变异函数,半变异函数在平面上的等值线是一组椭圆(图1-37)。椭圆的短轴和长轴称为主方向Principal Directions。对应于槛值的等值线上的每一点r到原点的距离是在o→r方向上半变异函数的变程。所以,对应于槛值的等值线椭圆称为各向异性椭圆,它是影响范围的一种表达。若平面为水平面,各向异性椭圆的长轴方向一般与矿体的走向重合(或非常接近)。因此即使矿体的产状是未知的,通过各向异性分析也可以确定矿体的走向。在三维空间,各向异性椭圆变为椭球体。 q 2区域各向异性区域各向异性的特点是半变异函数的槛值与变程均随方向变化,如图1-38所示。 q o o h a1 a2 C0 C0C p-p方向 q-q方向 p p r 图1-37 几何各向异性示意图 x h o C 图1-38 区域异性示意图 方向1 方向2 六、半变异函数平均值的计算 应用地质统计学方法进行参数估值时,需要计算半变异函数在两个几何体之间或在一个几何体内的平均值。设在区域Ω中有两个几何体V和W,如果在V中任取一点z,在W中任取一点z ,,z与z ,之间的距离为h,那么半变异函数在两点上的值为γh,记为γz,z ,。半变异函数在V和W之间的平均值就是当z取V中所有点、z ,取W中所有点时, γz,z ,的平均值,即 1-65 上式积分可以用数值方法计算。将V划分为n个大小相等的子体,每个子体的中心位于zii1,2,,n;同理,将W划分为n ,个子体,每个子体的中心位于z ,jj1,2,,n ,。这样,上面的积分可用下式逼近 1-66 当V和W是同一几何体时,即为半变异函数在几何体V内的平均值 1-67 式中,zi和zj都是V中的子体中心位置。 根据半变异函数的定义,γzi,zj,,xi 和xj,分别为区域化变量X在zi 和zj,处的取值。这样,式(13-66和1-67也可分别改写为以下的形式 1-68 1-69 如果几何体W代表的是一个样品,用ω表示,样品的中心位于z0,样品值为x0,而且样品ω的体积很小,不再划分为子体,即n ,1,那么式1-66变为 1-70 式1-68变为 1-71 称为半变异函数在样品ω与几何体V之间的平均值。 如果V也代表一个样品ω,,ω,的中心位于z0,,ω,的取值为x0,,ω,的体积很小,不再划分为子体n 1,那么式1-70和1-71分别变为 1-72 1-73 称为半变异函数在两个样品之间的“平均值”。 当然,如果样品的体积较大,需要把样品也划分为子体时,半变异函数在样品与几何体之间,样品与样品之间的平均值和半变异函数在两个几何体之间的平均值是一回事。因此,式1-66或1-68是计算半变异函数平均值的一般公式,其它公式都是这两个公式的特例。值得注意的是,当把样品也划分为离散点进行计算时,意味着半变异函数不是由具有一定体积的样品数值得来的,而是从无限小的点值得到的。这样的半变异函数称为点半变异函数Point semivariogram。但在实践中点半变异函数是未知的,半变异函数是通过具有一定体积的样品数据建立的。因此,在实际计算中一般把样品看作是“不可再分”的。 七、克里金Kriging 由于地质统计学法的基本思想是由Danie Krige(丹尼克里金)提出的,所以应用地质统计学进行参数估值的方法被命名为克里金法Kriging。克里金估值是在一定条件下具有无偏性和最佳性的线性估值。所谓无偏性就是参数估值与真值之间的偏差的数学期望为零,即 1-74 所谓最佳性是指估计值与真值之间偏差的平方的数学期望达到最小,即 1-75 也称为估计方差Estimation Variance,用表示;用克里金法进行估值的估值方差称为克里金方差Kriging Variance或克里金误差Kriging Error,用表示。 所谓线性估值是指未知量的估计是若干个已知取样值xi的线性组合,即 1-76 式中,bi为常数。 设从区域Ω中取样n个,样品ωi的值为xii1,2,,n;Ω中的一个单元体V的未知真值为(图1-39。那么,用这n个样品对的克里金估值即为式1-76。 1 ω V x1 x2 ω ω 2 xi i 图1-39 克里金法示意图 根据无偏性要求,有; 如果在区域Ω内,区域化变量满足内蕴假设且“无漂移”, E[xi]E[]μ(μ为参数在Ω的平均值),上式变为 ,消去μ得 或 1-77 因此,估值具有无偏性的充要条件是取样值的权值之和为1。 将几何体V看作是由m个相同的子体组成,每个子体的值为,那么等于子体值的平均值,即 。这样克里金方差为 令xj’,表示另一个样品ω,j的值(当j ’i时,xj’,xi);vl 表示V中的另一个子体的值(当1=k时,v1 vk。那么当时上式可以改写为 1-78 从公式1-73可知,即半变异函数在样品和间的平均值;从公式1-69可知,即 半变异函数在几何体V内的平均值;从公式1-71可知,即半变异函数在样品ωi和V之间的平均值。 将这些等价关系代入式1-78,得 1-79 这样,最佳估值就是在的条件下求达到最小值时的权值bii1,2,,n。应用拉格朗日乘子法,得拉格朗日函数 1-80 式中,为拉格朗日乘子。在的条件下求达到最小的条件是拉格朗日函数对bii1,2,,n和λ的一阶偏微分为零,即 对 1-81 将式1-79代入式1-80中并求导,得1-81的具体形式为 1-82 将1-82展开,得 1-83 式1-83是有n1个方程组成的线性方程组,称为克里金方程组。解这个方程组就可求出n1个未知数b1,b2,,bn,λ。将求得的b1,b2,,bn代入式1-76即得的无偏、最佳估值。式1-83也可用矩阵形式表示,令 式1-83可改写为 AB=C1-84 上式的解为 B=A-1C1-85 如果区域化变量在Ω中满足二阶稳定性假设,。因此,克里金方程组也可用协变异函数表示,这时 式中,为协变异函数在两样品间的平均值;为协变异函数在样品ωi和几何体V之间的平均值。和的计算公式与和的计算公式在形式上相同。 [例1-4]如图1-40所示,矿床Ω中有一个方块,其边长为3,一个样品ω1 位于方块的中心,其品位为x11.2;另一个样品ω2 位于方块的一角,其品位为x20.5。半变异函数的方程为 假设品位在矿床中满足内蕴假设且无漂移,试用克里金法估计方块V的品位。 3 1 2 3 ω1 1.0 5 4 3 6 0.5 8 7 9 ω2 V h 1 2 图1-40 例1-4中方块与样品相对位置及半变异函数示意图 解 1计算 ω1与ω2之间的距离为,ω1和ω2到自身的距离为 0根据公式13-72,有 2计算和 将V等分为如图所示的9个小块,根据公式1-70,有 式中z0为ω1的位置,zi为第i个小块的中心位置。当i2时,样品ω1距第二小块的距离为1,因此γz0,z2γ10.5。类似地可以求出任意i的γz0,zi 仿照上述做法,求得 3建立克里金方程组并求解 将上面求得的和代入式1-83,有 解上面的方程组得 b10.673, b20.327, λ0.209 4求方块V的品位 方块V的品位的估值为 5计算方块品位的克里金方差 计算克里金方差需要计算,根据式1-67 式中zi和zj为两个小方块的位置,对于每一对zi和zj,γzi , zj 的算法与第2步中γz0 , zi 的算法相同。这里不再重复,只给出计算结果如下 0.683 将有关数值代入式1-79得 0.175 对于图1-13所示的样品分布及品位值,通过不同方向的实验半变异函数计算和拟合,得知半变异函数在台阶上是各向异性的,长轴主方向为135o(即东南西北),短轴主方向为45o(即西南东北),两个主方向上的半变异函数均为球函数 长轴方向N0.015 C0.200 a 137米 短轴方向N0.010 C0.160 a122米 用克里金法对台阶上所有方块的估值结果如图1-41所示 第九节 影响范围 当对块状模型中每一单元块的品位进行估值时,无论用什么方法,均需要确定由哪些样品参与估值运算。一般地讲,对被估单元块有影响的样品(即落入影响范围内的样品〕都应参与估值运算。 地质统计学把品位看作是区域化变量,而且用半变异函数描述品位在矿床中的相互关联特征。因此地质统计学为帮助确定合理影响范围提供了理论依据。如前所述,在大多数情况下品位的半变异函数的数学模型为球模型。球模型的特点是半变异函数γh随距离h的增加而增加,当h增加到变程a 时γh达到槛值。由于槛值为样品的方差,这表明当h≥a 时样品值变为完全随机的,样品之间失去了相互影响。因此,半变异函数的变程a 可以被看作是影响半径的一种表述。然而当存在各向异性时,不同方向上的半变异函数具有不同的变程,影响范围是一椭球体,即各向异性椭球体。因此,确定合理影响范围首先要建立各个方向的半变异函数,进行各向异性分析。 影响范围在品位、矿量计算中起着相当重要的作用,在某些情况下,所选取的影响范围不同,矿量计算结果会有很大的差别。然而,确定合理影响范围是一件不容易的事,需要对矿床的成矿特征有深入的了解,同时也需要丰富的实践经验。地质统计学可以帮助确定合理的影响范围,但并不意味着各向异性椭球体(在各向同性情况下为球体)就是最合理的影响范围,最后决策应是综合考虑各种因素(包括经验)的结果。 600N 0 0 0 0 18 18 18 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 18 18 18 42 42 49 49 49 0 0 0 0 0 0 0 0 0 0 0 0 18 18 30 42 42 49 49 49 0 0 0 0 0 0 0 0 21 21 21 40 40 52 58 53 39 42 45 45 0 14 14 14 0 0 0 0 21 21 21 40 40 40 69 69 38 38 43 43 0 14 14 14 0 0 0 0 21 30 32 37 35 48 70 67 63 67 70 65 48 35 14 14 0 0 0 0 0 39 39 39 32 32 72 72 81 81 89 89 48 48 0 0 0 0 0 0 0 39 39 31 27 37 77 73 69 66 79 78 79 83 89 63 9 9 9 9 0 0 10 16 17 23 83 68 83 45 72 72 101 101 89 89 9 9 9 9 300N 0 0 10 14 10 13 33 44 99 83 83 91 120 109 64 48 8 7 6 6 0 0 10 33 43 44 2 2 136 138 92 92 134 134 52 52 7 7 4 4 0 0 0 64 64 64 2 45 62 59 73 93 152 132 68 61 34 28 4 4 0 0 48 48 50 49 41 28 28 28 64 64 161 161 76 76 48 48 3 3 0 0 48 48 45 41 41 21 20 20 25 35 80 78 92 73 18 13 3 3 0 0 48 35 31 29 41 17 17 17 6 6 40 40 100 100 1 1 0 0 0 0 0 23 23 23 0 20 20 19 15 15 15 24 59 50 1 1 0 0 0 0 0 23 23 23 0 22 22 22 19 19 3 3 39 39 0 0 0 0 0 0 0 0 0 0 0 22 22 22 21 19 3 19 39 39 0 0 0 0 0 0 0 0 0 0 0 0 0 23 23 23 0 0 0 0 0 0 0 0 300E 600E 0 0 0 0 0 0 0 0 0 23 23 23 0 0 0 0 0 0 0 0 0 图1-41 基于图1-13数据的克里金估值结果 小 结 本章介绍了取样数据的处理与分析以及品位(矿量)估算的几种常用方法。地质统计学法被大部分人认为是一种较好的品位估值方法,尤其适用于品位变化大,矿岩界线由品位控制的矿床,如侵染型的贵金属矿床。必须指出,无论方法在数学上如何复杂,在理论上如何先进,方法本身不能增加已知信息量,而品位与矿量估算的准确程度主要取决于已知信息量的大小。另外,不同的估值方法有不同的适用条件。对于一给定矿床,方法A最为合适;而对于另一矿床,也许方法B最为合适。各种方法在不同矿床的实际应用经验是选择估值方法的重要依据。 在实践中,国际上通常的做法是在对矿床进行初步评价或是数据量不足时,选用较简单的方法(一般是距离N次方反比法)。当有了足够的数据,对矿床进行正式可行性评价时,选用地质统计学法。地质统计学法也常被用于开采过程中局部范围内的矿岩重新圈定和品位、矿量计算。 49