Aligned Indicator Conditional Simulation of Probability of Karst-Fissure Media in Karst Area of Northern China
-
摘要: 针对当前岩溶学研究领域中空隙位置定量预测和数学描述的难题, 选择研究程度高、资料丰富、具有典型代表性的山东济南泉域的含水岩溶裂隙介质, 运用序列指示条件模拟方法研究空隙间的连通性能, 预测岩溶裂隙密集带和强径流带的位置.研究表明, 它既可给出其模拟数值的波动区间, 且不同分位数值模拟图均能明显一致地指示出岩溶裂隙密集带和强径流带位置, 尤其是当样本点较少时, 也可获得几乎与样本点较多时的一致成果, 从而为该方法的推广应用提供了广阔的领域.Abstract: Aimed at the current difficult problem of fix quantitative forecast and math description about pore position in karstology, the water-bearing karst-fissure medium in Jinan spring area was selected as a typical research area because of its higherinvestigation degree, abundant data and typical representation. The connectivity among pore was studied and the position of karst-fissure concentration and strong runoff zones were forecasted by applying sequence-indicator conditional simulation. It shows that the method can give the fluctuation interval of simulated data, and at the same time, the obvious conformance of different quantile simulation map indicates the position of karst-fissure concentration and strong runoff zones, especially when sample size is small(the identical effect can be obtained as large sample size). Thereby, the method can be applied to a wider extent.
-
岩溶形态与其他地质现象(变量) 一样, 它的形成、产生受其所处的地质背景条件决定, 具有必然性和遗传继承性, 即具结构性.同时也有受内外界条件作用而激发突变产生的随机性[1~4].在地质统计学中, 为了刻划各类地质现象(变量) 具有的结构性与随机性的二重性, 引入区域化变量的概念.区域化变量是在采用随机函数来描述地质变量的结构性和随机性的基础上提出的[5].近30年来发展起来的地质统计学是研究区域化变量的有效工具, 特别是各种模拟方法可以很好地刻划地质变量在全区范围内的连续性和变异性[6~8].
序列指示模拟, 地质统计学模拟在国外已有深入的研究和应用, 它可以模拟出多种结果, 充分反映地质变量的空间连续性、变异性以及整体范围的不确定性.目前, 国内外地质统计学方法研究表明, 地质统计学模拟在水文地质和岩溶学领域有广泛的应用前景, 但各项研究均处于尝试阶段.如何发挥各种统计学方法的特点, 如何应用定性地质条件信息和定量信息进行空隙空间位置、特征参数值的估计、模拟, 尚缺乏系统性的方法研究.
本文运用序列指示模拟方法对济南含水岩溶裂隙系统结构的连通性进行定量模拟.
1. 序列指示模拟
设区域化变量Z (x) 满足二阶平稳假设, 即数学期望存在, 且等于常数E[Z (x)]= m; 协方差函数或变差函数只依赖于滞后h, 而与点的位置无关, 则条件模拟计算公式为:
式中: Zsc为条件模拟值; Zs为非条件模拟值; Zk* (x) 为真实值的克里格估值; Zsk* (x) 为非条件模拟值的克里格估值.
由于Zs (x) 与Z (x) 同构, 求克里格估值时数据构形又一样, 故两者克里格方程也一样, 解也一样, 所以有相同的权系数λα (α=1, 2, …, N), 依据
可将条件模拟计算公式改写为:
条件模拟的主要优点是能重复再现区域化变量真实值的离散性、波动性大小.在未观测点上也可以给出模拟值, 且只要重复进行, 就可得到许多个模拟的实现.但必须牢记: (1) 模拟值并非真实值; (2) 模拟值也不是最优的估计值.以Zsc (x) 作为Z (x) 的估计量时, 其估计方差是以Zk* (x) 作估计量时估计方差(即克里格方差) 的二倍.
此外, 模拟还可通过利用模拟点附近的样品集合或其他已知信息, 建立模拟位置x点处的不确定性模型——概率分布函数的方法来进行.建立不确定性模型是对Z (x) 的空间概率分布构造另一个相同概率模型的统计过程[9].
序列指示模拟是Journel①②近年来在指示克里金基础上提出的一种新的模拟方法.它与指示克里金一样可将样本点特异值对空间结构的影响给予削缓, 而又保持其随机性分布的特征.序列指示模拟的关键是要通过指示克里金构造每一个模拟(待估) 点的概率分布函数.序列指示模拟技术的基本步骤如下:
① Journel A G. Fundamentals of geostatistics in five lessons. 1994.
② Journel A G. Focusing on spatial connectivity of extreme-values attributes: stochastic indicator models of reservoir heterogeneity' s, SPE.1996.
(1) 设区域化变量Z (x) 有n个样本值{Z (xα), α=1, 2, …, n}, 记
在Zmin和Zmax之间取L个阈值且满足Zmin=z0 < z1 < z2 < … < zL < Zmax据此可将样本数据分成L+1组.一般L=3~9, 阈值的取法应保证各组所含样品的数量基本均匀.
(2) 定义z的阶梯函数(亦称指示函数)
对每个一般数据Z (xα), 用指示函数进行编码, 构成一个由0和1组成的指示数据列.这个数据列的特点是随着阈值的增加单调不减.
(3) 对于固定的每个阈值zl计算指示实验变差函数
再用通常的方法拟合适当的理论变差函数模型, 进而还可计算协方差函数.
(4) 对于N个待模拟的点, 随机地确定一个模拟顺序.
(5) 按选定的顺序, 对每一个模拟节点xj, 做L次指示克里金运算, 以确定Z (xj) 的条件累计分布函数值①②.
这样计算出的分布函数应该是单调不减的.若计算出的结果不满足这个条件, 需做适当调整.若
(6) 利用蒙特卡罗法获取xj处模拟的现实Zs (xj), 即从[0, 1]间的均匀分布中求取随机数P (P亦为累积概率值), 由概率分布函数F (xj, z) 的逆函数可求出相应P分位数, 即xj处的一个模拟值Zs (xj) =F-1 (P, xj).
(7) 将所得的模拟值Zs (xj) 加到已知值中, 即原来的已知点增加一个点.回转到第5步按第4步选定的顺序, 用同样的方法计算下一个点的模拟值Zs (xj+1).直至所有的模拟点计算完毕.如果需要再进行一次模拟, 可从第4步重新开始计算.
应说明的几点是: (1) 在每一个模拟位置估计概率分布函数时仅仅依靠该估计点附近的数据点.随着估计点位置的改变数据点也改变, 且有一些邻近已获得模拟值的新点值加入到数据集中.这意味着对于每一个模拟位置都需要解一个新的指示克里金方法集, 而且计算量会逐渐增加.在实际上, 仅用邻近的若干个点计算新的模拟值. (2) 因风险分析和前景估量需要很多模拟结果, 如果需要可进行多次模拟计算.
2. 研究区含水岩溶裂隙系统特征
济南泉域属华北型地层区, 南部山区自南而北分布: 太古界泰山群(AR) 变质岩系, 古生界寒武系(∈) 灰岩及页岩, 中、下奥陶统(O) 灰岩、泥质灰岩及白云质灰岩, 地层产状走向近东西, 倾向北, 倾角平缓.其中, 中、下奥陶统至山前逐渐隐伏于第四系(Q) 之下.市区北部地下分布中生界燕山期(ν) 辉长岩、闪长岩; 东郊和西郊的北部地区下伏古生界石炭系(C)、二叠系(P) 砂页岩与中奥陶统灰岩成假整合接触.新生界第四系松散堆积层广泛分布于山前倾斜平原, 其中西郊分布有玉符河和北沙河形成的冲积-洪积层.区内在地质构造上为一单斜构造, 断裂有三组: 北北西向的东坞断裂、千佛山断裂、石马断裂、平安店断裂和马山断裂; 北北东向的炒米店断裂; 北东向的港沟断裂, 这些断裂将济南单斜分割成三大断块.即东梧断裂-千佛山断裂断块、千佛山断裂-平安店断裂和平安店断裂-马山断裂断块(图 1).
济南泉域位于北方半干旱岩溶区, 岩溶形态以溶隙和溶孔组成的岩溶裂隙脉状网络系统发育为特征, 其发育方向主要受构造裂隙系统控制.因此, 选取济南泉域岩溶裂隙作为反映该研究区岩溶发育特征的区域化变量.在对泉域岩溶裂隙介质特征分析的基础上, 根据岩溶裂隙点的野外测量的构造裂隙参数, 按正态分布的“一倍方差”原则, 将极宽裂隙和个别散乱点作为风化、卸荷裂隙等随机因素滤掉, 参与计算的构造裂隙合格率达到85%以上, 共计122个岩溶裂隙实测点.
将经处理后的裂隙、岩溶参数代入上述有关公式, 计算出各点各向同性化的岩溶裂隙渗透系数.并在充分分析该区岩溶裂隙发育特征的基础上, 选取渗透系数作为反映含水岩溶裂隙系统的区域化变量.
3. 岩溶裂隙渗透系数的序列指示条件模拟
济南泉域岩溶裂隙渗透系数的序列指示条件模拟均值等值线图(图 2) 显示了岩溶构造场对渗流场的控制作用.即区域上渗透系数等值线呈近东西向展布, 从南向北由补给区到径流排泄区渗透系数逐渐增高, 特别在山前地带有一近东西向的高渗透系数带.
图 2 渗透系数序列指示条件模拟均值等值线(122个样本点, 渗透系数单位m/d)Fig. 2. Sequence-indicator conditional simulation of permeability coefficient of karst-fissure medium in studied area序列指示条件模拟分位数等值线图(122个样本点) (图 3) 表明, 济南泉域含水岩溶裂隙系统具有很强的非平稳性.即分位数为5%, 25%, 50%和75%对应的模拟均值等值线图所反映的基本形态都表明中部存在一条南北延伸、连通性极好的岩溶通道.西北角虽有一高度最大的溶洞, 但其基本处于孤立状态.东北部溶洞向西南逐渐归并于中部岩溶通道.虽然按不同的分位数可给出不同的模拟值, 但反映总体连通性的规律不变.
将122个已知点与88个已知点渗透系数的序列指示条件模拟成果图(图 2、图 4) 进行对比分析, 其模拟结果基本一致, 说明样本点稍少时序列指示条件模拟方法也可取得较理想模拟效果.不同分位数等值线图虽然模拟出的渗透系数大小有变化, 但它们分布的形态、高渗透系数值的位置连通方向均趋于一致.表明序列指示条件模拟方法是定量模拟岩溶空隙连通性的好方法.
图 4 岩溶裂隙介质渗透系数序列指示条件模拟E型值等值线(88个样本点, 渗透系数单位m/d)Fig. 4. Sequence-indicator conditional simulation of permeability coefficient of karste-fissure media in research area图 2~4都突出反映了北北西向和山前东西向渗透系数高值带条的出现, 这不仅刻画出该泉域地下水补给、径流、排泄区的分布, 地下水的总体流向及主径流带的位置, 而且指示出构造场对渗透场的控制作用.沿北北西向千佛山断裂中、北端发育一条主要的岩溶裂隙密集带, 它也是地下水的主径流带.北部山前同样存在一条近东西向展开的岩溶裂隙密集带, 且处于排泄区内, 渗透系数大, 地下水相对丰富.
4. 结论
针对当前岩溶学研究领域中空隙位置定量预测和数学描述的难题, 本文选择研究程度高、资料丰富、具有典型代表性的山东济南泉域的溶隙介质, 应用序列指示条件模拟方法研究空隙间的连通性能, 预测岩溶裂隙密集带和强径流带的位置. (1) 序列指示模拟是从随机角度出发探讨问题, 认为模拟点均遵循面上样本点的概率分布函数.序列指示条件模拟通过多重指示克里格方法考虑周围已知点和已模拟出的点对该点的影响, 即区域结构成分.再根据模拟点的概率分布函数, 求取不同概率条件下的随机值.因此序列指示条件模拟反映模拟点上渗透系数波动范围, 从整体角度则显示岩溶裂隙介质渗透性能. (2) 序列指示条件模拟是模拟岩溶裂隙介质连通性能, 指示或预测岩溶裂隙密集带、强径流带位置的行之有效方法.它既可给出其模拟数值的波动区间, 且不同分位数值模拟图均能明显一致地指示出相对高值位置, 尤其是当样本点较少时, 也可获得几乎与样本点较多时的一致成果, 更为该方法的推广应用提供了广阔的领域.
-
图 2 渗透系数序列指示条件模拟均值等值线(122个样本点, 渗透系数单位m/d)
Fig. 2. Sequence-indicator conditional simulation of permeability coefficient of karst-fissure medium in studied area
-
[1] 袁道先. 中国岩溶学[M]. 北京: 地质出版社, 1993.YUAN D X. Karst of China[M]. Beijing: Geological Publishing House, 1993. [2] 卢耀如. 中国岩溶——景观·类型·规律[M]. 北京: 地质出版社, 1986.LU Y R. Karst in China—landscape·types·rules[M]. Beijing: Geological Publishing House, 1986. [3] 杨银湖, 黄正发. 京珠线湖北省南段岩溶地质问题与勘察对策[J]. 地球科学——中国地质大学学报, 2001, 26 (4): 361-364. https://www.cnki.com.cn/Article/CJFDTOTAL-DQKX200104007.htmYANG Y H, HUANG Z F. Karst geological problems and their countermeasures in southern Hubei section of BeijingZhuhai speedway[J]. Earth Science—Journal of China University of Geosciences, 2001, 26(4): 361-364. https://www.cnki.com.cn/Article/CJFDTOTAL-DQKX200104007.htm [4] Chin-Fu Tsang. 非均质介质中地下水流动与溶质运移模拟——问题与挑战[J]. 地球科学——中国地质大学学报, 2000, 25(5): 443-450. https://www.cnki.com.cn/Article/CJFDTOTAL-DQKX200005000.htmChin-Fu Tsang. Modeling groundwater flow and mass transport in heterogeneous media: issues and challenges[J]. Earth Science—Journal of China University of Geosciences, 2000, 25(5): 443-450. https://www.cnki.com.cn/Article/CJFDTOTAL-DQKX200005000.htm [5] 王仁铎. 线性地质统计学[M]. 北京: 地质出版社, 1993.WANG R D. Linear geostatistics[M]. Beijing: Geological Publishing House, 1993. [6] 潘国成. 地质统计学中的区域化变量理论[J]. 世界地质, 1997, 16(2): 85-93. https://www.cnki.com.cn/Article/CJFDTOTAL-SJDZ702.016.htmPAN G C. Regionalized variable theory for geostatistics[J]. World Geology, 1997, 16(2): 85-93. https://www.cnki.com.cn/Article/CJFDTOTAL-SJDZ702.016.htm [7] 肖斌, 赵鹏大, 周龙茂. 归来庄金矿床w(Au)/w(Ag)异常的地质统计学研究[J]. 地球科学——中国地质大学学报, 2000, 25(1): 79-82. https://www.cnki.com.cn/Article/CJFDTOTAL-DQKX200001016.htmXIAO B, ZHAO P D, ZHOU L M. Geostatistical study of w(Au)/w(Ag)anomaly in Guilaizhuang gold deposit[J]. Earth Science—Journal of China University of Geosciences, 2000, 25(1): 79-82. https://www.cnki.com.cn/Article/CJFDTOTAL-DQKX200001016.htm [8] 潘国成. 地质统计学中结构分析的理论与方法[J]. 世界地质, 1997, 16(3): 70-82. https://www.cnki.com.cn/Article/CJFDTOTAL-SJDZ703.013.htmPAN G C. Theory and technique of structural analysis for geostatistics[J]. World Geology, 1997, 16(3): 70-82. https://www.cnki.com.cn/Article/CJFDTOTAL-SJDZ703.013.htm [9] 潘国成. 地质统计学中的估值技术和条件模拟[J]. 世界地质, 1997, 16(3): 83-100. https://www.cnki.com.cn/Article/CJFDTOTAL-SJDZ703.014.htm -
下载:
下载:
点击查看大图
计量
- 文章访问数: 3471
- HTML全文浏览量: 116
- PDF下载量: 5
- 被引次数: 0
