petrel学习-随机模拟原理部分

第一节类特征的未知数据作一种线性无偏、最小方差估计(BestLinearUnbiasedEstimator(BLUE))一、克里格的基本思想及数学表达V(

V(x0)VV

Zx)dx(简记为

ixivi上的的平均值Zv(xi(简记为Z)niivi(i1,2,nV1-1V V

图1- n iZ* iiVVV计方差最小的前提下n个权系数i在这样的条件下求得的i所构成的估计量Z*V E[Z(xm(常数)Y(xZ(x)mE[Y(x0E[Y(xYyC(xyY

Y(x)dxV

Z(x)dxm

Y*Y，其中YZm(i1,2, ii 故只要求得Y的估计值Y*，就能得到Z的估Z*。Y*是Y

E(Y*)E(Y)E(Zm) i

iE(Y*E(Y，即Y*是Y 为了求出使估计方差2E[YY*2为最小时的权系数(i1,2,n 2E[YY* E[Y2]2E[YY*]E[Y*2 E[Y(x)Y( 2

E[Y(x)Y(V2 iV V iiE[Y(xi)Y(xji1j

i C(x,y)dxdy

C(x,x)dx

C(x,xV2

iV

i V 2C(V,V)

i i1j C(x,V)C(x,x i i1i

EE

n2C(xi,V)2jC(xi,xj) (ijnjC(xi,xj)C(xi,V (ij

从这个方程组中解出i(i1,2,n，即为所求的简单克里格权系数，它必定满足最 ijC(xi,xj)iC(xi,V i1j i2C(V,V)C(x,V iEEKKi(i1,2,n从方程(1.18)中解出i之后，即得到YV的简单克里格估计Y* jjZ*mY*mYm

(Z jj

j K Z*jZm1j Kj j 普通克里格方程是在区域化变量Z(x)的数学期望未知的情况下，建立的克里格方程，故此时E[Z(x)]m为未知常数，若要求估计量是无偏估计，就必须有限制条件，称为无ni1n

C(x,x) C(x,V),(i1,2,,jjj

ijC(xi,xj)iC(xi,V)i1j i 2C(V,V)C(x,V)C(x,x

i1

C(V,V)C(x,V) iK2C(V,V)[]T[M]K))

(x,

)(x,V),(i1,2,,jj jii

2(x,V)(V,V) i的平均结果，则式(1.23)中样点之间的协方差C(xi,xj)C(vi,vj)，这积vi就称为样点xi的支撑或承载(support)，如果对支撑进行克里格估计j

j(vi,v

)(vi,V),(i1,2,,

ii (v,V)(V,V) i只有当协方差矩阵[C(xi,xj)]nn(即矩阵[K]的左上角nn阶方阵)是严格正定的，的点协方差函数C(h是正定的(若用变差函数(h，则要求(h是条件正定的)，且数据合，则由克里格方程组给出的Z Z(v)，及 对于(1.23)和(1.25)所用到的协方差函数C(h)和变差函数(h)的模型，无论它们所表普通克里格矩阵[K，只取决于样品支撑的几何特征(空间位置)，而完全不依赖于待估程组的系数矩阵[K]也相同。从而只需求一次逆矩阵[K]1。若估计构形(待估支撑与全体样品支撑的构形)也相同，则矩阵[M]也不变。即只需解一次克里格方程组，就可以得到线性估计量的权系数i(i1,2,n。这就大大的节省了计算时间。考虑到这一点，在实际数据构形的几何特征(xi,xj信息样品支撑viV之间的距离(vi,V反应区域化变量Z(x)空间结构特征的变差函数模型二、各种克里格方个残差的和，所以“具有趋势模型的克里格(KT)Z(u)分成两部分∶Z(u)m(u) 其中m(uE[Z(u，是趋势部分，通常是用一个光滑的确定性函数来模拟，或用拟Km(u)akfkk上为具有限制条件的正规方程组，称为KT线性系统。KT的估计值为∶nZ*(u)(KT)(u)Z(u

(KT)(u)C(uu)K(u)f(u)C(uu1 k 1,2,, (KT)(u)f(u)f k0,1,, 1 其中(KT)是KT权值 fk(u，例如∶如果认为对z(u)的空间起作用的是一个周期性的成分，那么应该考虑具有某种周期和相位的正弦函在无法得到有关趋势的形状信息时，如何根据式(1.33)z数据两分成趋势部分和残差解释的情况下，通常选用u的低阶(2)多项式来表示趋势，其中u(x,y∶m(ua045m(ua0a1xm(uaaxayax2ay2a z数据时，应根据过滤掉趋势m(uz数据的线性组合来推断残差的协方差函数CR(hz(uhz(u0阶的趋势m(u)a0z(u2h2z(uhz(u的二阶差分将过滤掉任何一阶趋势m(u)a0a1u，等等m(ua0a1f1uf1u项与一个二级(外部)变量相等，二级变量的光滑性被认为与所要估计的主变量Z(u)有关。E[Z(u)]m(u)a0a1 y(u)反映了z变量的空间趋势(对应于两个参数a0和a1)且f1(u)y(u)，即∶nz*(u)(KT)(u)Z(u n(KT)(u)C u)(u)(u)y(u)C(uu1 1,, (u)y(u)1 其中 是克里格(KT)权值，是拉格朗日参数关系(1.36)必须具有一定的物理意义。例如∶如果二级变量y(u)代表的是到一个反射面z(u)应该(在平均意义上)与旅行时成比例，因此关系式(1.36)就有意系统(1.37)中，变量Z(u)和Y(u)之间的互协方差函数不起作用，这与协同克里格不同。Z(u)Z0(u)Z1(u)ZLLCZ(h)Cl

例如∶(L+1)RFZl(u)z样品的(L+1)个套合协方差函数来模拟。 Zl da(u)Z(u nd(u)C(uu)(u) C(uu

,n 1d(u)

(1.41)E{Zl(u0,l0,l01的条件下估计值(1.40)根据孔隙样品数据和相关的声波数据估计孔隙度值z(u)。 (u)

(u)Z(u

2)2

2 其中是作用于nz是作用于ny 模型，它包括Z协方差函数CZ(h)，Y协方差函数CY(h)，互Z-Y协方差函数CZY(h)Cov{Z(u),Y(uh)}Y-Z协方差函数CYZ(h)z未知量相关性较好的数据(往11 (u)11

(u)Z(u) CZY(h)，且可以使用下面的CZY(h)BCZ(h), 其中B CY(0)/CZ(0)ZY(0)，CZ(0),CY(0)是Z和Y的方差值ZY(0)，是同位小，且要记住y只不过是二级变量。z(u)y*ulnz(u，然后利用简单克里)z(u)的估计值；尤其是数变换ey*u)Z(u)的有偏估计,简单对数克里格的估计值 2(u)z(u)expy(u) 法，尤其是MG和IK方法，在随机模拟中得到了发展。MG方法中，把SKccdf均值和方差并没有反被ccdf中提取的正态得分的模拟值反转换回去，这种转换是不需要类似地，IK得到的ccdf估计值并没有被反转换回去，而是直接从ccdfz变函数，对正态得分数据进行SK或OK。第二节随机建模的原理和方法的、确定的特征和性质（Haldorsen和Damisieth,1990），同时又具有复杂性。这同。既有在一个小岩心栓上的直接测量数据，也有勘探低分辨率的间接的测工反映Z(u)空间分布的可供选择的、等概率高分辨率实现。实现｛Z(i)(u)，1、两步法模拟模型Z(u)不能模拟严重非均质性，如穿越岩石类型之类的物理边界。2、多变量联合模拟相关性最好的变量，Journel等称之为原始变量（主要变量），然后用条件分布3、序贯模拟方法许多模型的近方法是在某一位置u处给定最相关的变量的条件下，从条件有点，这些点包括原始的样品点和已模拟好的点。假设N个随量Zi的联合随量Zi可以指研究区内N个节点样本的参数，或在同一位置处测得的N个不同的变量，或代表N个节点K个变量，式中“|”表示“在某一条件下”，条件分布表明，要从累积条件分布函数取出N个变量的样本来，可以用N个连Z1的累积条件分布函数上取一个值Z1(1)，该值是在给定的n个原始样品数据条的n个增加到n+1个。同样在给定的(n+1)个数据的条件下，再从Z2的单变量的ccdf中，取Z2(1)1，从而变为(n+2)。就这样，一步一步地按顺序地考虑所有N个随量Zi，那么就可得到集合｛ZI(1),i=1,…,N｝，代表了N个相互关联的随量Zi的一个实现。如果还需要另一个模型，则需要再重复这样顺序。在这种序贯模拟过程中，需要确定NP｛Z1≤z1|(n)｝…P｛ZN≤zN|(n+N-序贯模拟方法是独立于建立的单变量ccdf方法和模型的，而在高斯方法中，所有的ccdf都假设为高斯分布，它们的均值和方差由N个简单的克里格系统来给出，而在序贯指示模拟中，ccdf是直接由指示克里格系统给出。数据也应加以考虑。用邻࿑的数据作为条件，在模拟中只能体现给定距离内的变量的统计特征。例如，条件范围必须大到足以体现变差函数。这就需要从1到Nn个节点的顺序最好是随机的。因为，如果n个节点是按行时，将会沿行出现人为的效应。4、误差模拟模拟的Z值将是唯一的估计值和模拟误差的迭合：ZR(u)必须是独立的，或至少要与估计值Z*正交。实际上，误差是独立地模拟后加到估计值如果估计Z*(u)是把Z(u)正交地投映到数据的某一空间中而得，在这种情况下，误差Z*(u)-Z(u)与估计Z*1是可以满足的。由于误差的协方差2Z(u)是平稳的，由于在某位置处u的数据特征在其它地方是不一样的，误差协方差不是平稳的，解决的办法某一非条件模拟的实现Z(1)(u)，在所有节点外和实际位置处都有模拟值，这一实现具有随机函数模型Z(u)同样的协方差。保留在n个原始数据点处的模拟值Z(1)(u),=1,…,n。再用估计方法对实际的数据进行计算。这就给出了估计值Z*(1)(u)的模拟。那么模拟误差为r(1)(u)=Z(1)(u)-Z*(1)(u)，模拟误差就可以简单地加到实际的估计值Z*(u)，也即：数据的，即Zc(1)(u)=Zc*(1)(u)，则对于任何的=1,…,n有Zc(1)(u)=Zc*(u)=Z(u)Z(1)(u)的变差函数与原始Z(u)的变差函数模5、P场模P场（概率场）基于原始数据的原的ccdf函数中提取Z实现。假设P(1)(u)和P(l)(u1)是u和u1区间均匀分布，协方差模型来自数据均一化的样品协方差。因为P场实现｛P(1)(u)，uA;L=1,…,L}是非条件模拟，故它们能由任何快速模拟算法来产生，如转向带、非条件分形或简单地随机窗口移动平均法设F(u;Z|(n))和F(u′;Z|(n))是u和u′处的ccdf函数，仅受限于n个原始数所以它有快速之优点。该ccdf函数可由Z连续数据正态得分转换实现的高斯克里格法得到，或由指示数据实现的指示克里格方法获폖。Z模拟值从这些函数用空间相关概率值P(1)(u)和P(1)(u)求取Z(l)(u)=F*—限制），但Z实现是条件的。在数据点u,ccdf函数F(uz|(n))具零变差，并集中于数据值Z(u)，因此不管概率值P(1)(u)如何，ccdf函数将总是返回数据值；Z(l)(u)=Z(u)。6、退火模拟模拟退火（SA技术是把最优化问题与统计力学热平衡问题进行类比得来的。性。几个不很严格的准则已在随机的逐次近似法(stochasticrelaxation)中得到应Boltzmannp(E)exp(E2E1这里T是某个热平衡系统的温度，EKb称为Boltzmann常数，波兹曼（Boltzman）T时，热平衡状态下，系统所具有的能量按概率分布于不同的能量状态中。Kb它是把温度和能量关联起来的常已出现的模拟退火方法有：Metropolis算法、热浴法、协同退火法、并行模拟退能量状态E1到能量状态E2变化的概率P(E)如果E2＜E1，则系统将变化，而建立初始的控制参数 和降低它的过程建立接受 交换1概率分P{accept}

e(OO)/

数Oold=Onew；控制参数t狄尤史(C.V.Deutsh)1992年提出了一种经验方法，该方法用7个参数来表达某一退火计划，并给出最佳的参数选择（1建议的模拟退火计划（ClaytonV.Deutsh*KS淬///23533表中：t0：初始的控Kmax:在一个控制参数下，进行扰动的最大次数（通常取节点数 倍），当达到Kmax时，用λK接爱：扰动被接受的次数，在每一控制参数下，当接受了K接受S：停止控制参数：如果达到KmaxS次，那么算法停止（通常S3）△0：指示目标函数收敛的一个参数，即一个小的目标函数值。（一）随机模型的基本分类机模拟的方法从不同的角度进行了分类，如HaldorsenDamseleth(1990)、AlabertModot(1992)、DeautchJournel(1992,1996)以及Srivastava(1994)。共（三）基于象元的随机模型1、高斯模拟方法 如果连续的空间随机函数

A是数目并不太多的一些具有类似由于高斯模拟中要求所模拟的随机函数的ccdf为正态分布，因此，整个模再将模拟结果反变换为区域化变量。对于符合正态分布的变换后随量，可以果Z积分布函数的Y数据；如果多变量高斯模型适用于Y（正态得分变换后样品数据），则ccdfY将模拟值Y(l)(u)序贯高斯模拟是应用较为广泛的连续变量的模拟方法。从算法来讲比较稳（如储层的渗透率2、截断高斯随机域F(u)=i如果ti-设ni中岩相，其指示值可用高斯随机函数来定义：I Yx 如果Yxi1,i 因此，点x属于第i种岩相当且仅当Yxi1,i。其中，i为截断值。nFxcodiIi1Yxi因 在位置n处取值I，当且仅当位置属于相I，即Ii1Yxi3、序贯指示模拟所谓指示变换就是将原始数据按照不同的门槛值，编码成1或0的过程，假设在X处的参数Z(x)，对门槛值为Zc的指示变换可写成：

I(u,z)

Z(u)1，0；也可以是地质家的解释和推断，因此，10它不是用来某个具体值，而是作指示。上式是一种等权指示法，对未知数通常也用不等权法: