CN104459776A - 一种断裂分形特征优化计算方法 - Google Patents

Abstract

本发明涉及油气田勘探开发领域，尤其是一种断裂分形特征优化计算方法。本发明在断裂带数据化后，通过设置边界参数、起始位置，定义统计箱子的参数，通过编写程序使统计箱子不断的移动，并且计算统计箱子所在位置的分维值，实现断裂分形特征的程序化计算。本发明对于断裂分布特征定量评价、断裂分形特征定量描述等多个方面具有较高的实用价值，并且预测成本低廉、可操作性强，大量减少人力的支出，评价结果在预测油气勘探开发重点区域、裂缝优势分布区域以及低级序断层的分布特征等多个方面有一定的指导意义。

Description

一种断裂分形特征优化计算方法

技术领域

本发明涉及油气田勘探开发领域，尤其是一种断裂分形特征优化计算方法。

背景技术

在断裂构造的研究中，人们已经认识到，断层、裂缝在几何学、运动学以及动力学特征上具有统计意义的自相似性，这是分形理论在构造地质研究中的体现。越来越多的研究表明，岩石破碎过程具有随机自相似性，断裂的分布、几何形态具有明显的分形特征。分维是定量评价断裂构造的一种准确、有效的指标，断裂的分维值不仅与断层的长度有关，而且与断裂的条数、长度以及平面组合特征等有关。断层和裂缝在成因上具有一致性，在相同的应力场背景下形成，因此确定二者之间的内在定量关系，依据断裂分布可以预测裂缝的发育规律。根据计算维数的定义，可分为容量维、信息维、相似维以及关联维等，在分形几何构造地质复杂性的应用研究中，容量维和信息维较实用。

断裂分形统计过程中，往往采用人工的手段，工作量大、费时费力，并且计算后结果人为因素较多，尤其在度断裂求取面积信息维、长度信息维时，断裂走向不规则、断裂宽度变化的随机性，这些因素都困扰裂缝分形特征准确-快速评价，因此，有必要展开利用计算机程序，描述断裂分形特征。

发明内容

本发明旨在解决上述问题，提供了一种断裂分形特征优化计算方法，它解决了断裂分形统计过程中，采用人工手段工作量大、费时费力且预测分维值不准确的问题。

本发明的技术方案为：一种断裂分形特征优化计算方法，具体步骤如下：

第一步根据地震解释构造图，对断裂带数据化

利用相关软件，对地震解释构造图中的断裂数据化，采取折线式数据化方法，如图2所致，分别在断裂的上升盘对应的断层线与下降盘对应的断层线依次取点，其中数据的起点与终点分别对应断裂的两端点，在断层走向变化大的区域，采用密集取点方法，保证取点次序的同时又要保证断裂的原始形态，使取点能最真实的反映断裂的形态，设定每条断裂读取的点序号标记为(a₀、a₁、a₂...a_n-1、a_n)。

第二步设置参数，包括起始位置、充填步长以及边界条件

根据研究区的范围，设置工区的起始范围、边界范围以及对断裂带充填步长b(单位m)。断裂带充填步长b主要与断裂带的宽度有关，一般可以设置为研究区小断裂宽度的十分之一，或者更小，一般断裂带充填步长b越小，计算断裂分维值精度越高，但计算机的工作量也越大。

第三步断裂带充填

如图2所示，每条断裂读取的点的序号依次标记(a₀、a₁、a₂...a_n-1、a_n)，对每条断裂，利用步骤一对断裂数据化后，每条断裂由一系列依次排列的三角形组成，对于任意的三角形，设其由三点a_m、a_m+1、a_m+2组成，其中m的变化范围0—(n-2)，点a_m、a_m+1、a_m+2对应的X轴、Y轴坐标分别为(x_m，y_m)、(x_m+1，y_m+1)、(x_m+2，y_m+2)，可以依据三个点坐标数值确定充填的范围-矩形ABCD(图3)，根据断裂带充填步长b，将矩形ABCD划分为一系列边长为b的小格子(图4)，依据断裂带充填步长b以及矩形ABCD的边界坐标值，可以表示矩形ABCD内任意一个小格子中心所在的坐标数值s(x_p，y_q)，通过公式(1)-(3)，可以判断任意的点s(x_p，y_q)是否落入断裂带(点a_m、a_m+1、a_m+2)组成的三角形内：

参数aa₁、bb₁、cc₁表示为：

$\left\{\begin{array}{c}{\mathrm{aa}}_1=(x_{m+2}-x_m)\left(y_{m+1}{-y}_m\right)-(y_{m+2}-y_m)(x_{m+1}-x_m)\\ {\mathrm{bb}}_1=(x_{m+1}-x_m)(y_{m+2}-y_m)-(y_{m+1}-y_m)(x_{m+2}-x_m)\\ {\mathrm{cc}}_1=(x_m-x_{m+2})(y_{m+2}-y_{m+1})-(y_m-y_{m+2})(x_{m+2}-x_{m+1})\end{array}\right.---\left(1\right)$

参数aa₂、bb₂、cc₂表示为：

$\left\{\begin{array}{c}{\mathrm{aa}}_2=(x_p-x_m)\left(y_{m+1}{-y}_m\right)-(y_{q+2}-y_m)(x_{m+1}-x_m)\\ {\mathrm{bb}}_2=(x_{p+1}-x_m)(y_{m+2}-y_m)-(y_{q+1}-y_m)(x_{m+2}-x_m)\\ {\mathrm{cc}}_2=(x_p-x_{m+2})(y_{m+2}-y_{m+1})-(y_q-y_{m+2})(x_{m+2}-x_{m+1})\end{array}\right.---\left(2\right)$

参数aa、bb、cc表示为：

$\left\{\begin{array}{c}\mathrm{aa}={\mathrm{aa}}_1\·{\mathrm{aa}}_2\\ \mathrm{bb}={\mathrm{bb}}_1\·{\mathrm{bb}}_2\\ \mathrm{cc}={\mathrm{cc}}_1\·{\mathrm{cc}}_2\end{array}\right.---\left(3\right)$

当aa≥0且bb≥0且cc>0时，认为小格子s(x_p，y_q)落在由点a_m、a_m+1、a_m+2组成的三角形内。依据公式(1)-(3)，在对断裂网格化的基础上，可以筛选出所有落入断裂带的小格子s，在此基础上，实现对断裂带的充填；落在断裂带内小格子s的总数与每个小格子的面积乘积就可以近似认为是断裂带的面积。

第四步统计箱子(F)参数设置

所述的统计箱子(F)就是指在断裂分形中所要计算-评价维数的网格单元；设置统计箱子(F)的边长L(正方形)，统计箱子(F)每次的移动步长为⊿l，同样移动步长⊿l越小，计算的研究区分维值精度越高，但计算机的工作量也越大；通过编写程序，使统计箱子(F)移动，这样可以统计边界范围内每个统计箱子(F)内断裂的分形特征。如图5所示，通过编写程序使统计箱子(F)首先沿X轴移动，每次沿X轴正方向移动步长⊿l，直至到达右边界，完成一次X轴方向移动；之后返回统计箱子(F)起始位置，沿Y轴正方向移动步长⊿l，开始第二次X轴方向移动，直至到达右边界，完成第二次X轴方向移动；之后返回统计箱子(F)起始位置，沿Y轴正方向移动2×⊿l；......直到统计箱子(F)沿Y轴达到上边界时，进行最后一次沿X轴方向移动，到达边界的东北边界，完成在边界范围内的移动。

第五步断裂面积信息维与容量维的求取

对于平面上的一个统计箱子(F)，将这个统计箱子(F)所在的矩形区域分割成若干个边长为r的网格，记F中的点落在小方格数目为N(r)，则定义容量维D_k为：

$D_k\left(F\right)=-\underset{r\→0}{\lim }\frac{\ln N\left(r\right)}{\ln r}---\left(4\right)$

信息维(D_I)：容量维只考虑了统计箱子F中小网格数目，却没有考虑不同边长的小网格内覆盖的点数，因此引入信息维D_I的概念：

$D_I\left(F\right)=-\underset{r\→0}{\lim }\frac{I\left(r\right)}{\ln r}---\left(5\right)$

其中，记P_i是统计箱子F中的点落在第i个小网格中的概率。

通过步骤三，筛选落入断裂带的所有小格子组成的点集T_sum，通过编写程序，统计落入每个统计箱子(F)中的数目、X坐标、Y坐标。然后将这个统计箱子(F)所在的矩形区域分割成若干个边长为r的正方形网格后，同样可以统计落入每个小网格的数目，通过不断变换r(一般设置为统计箱子边长L的1/3-1/10)不断统计落入每个小网格的数目。

当不考虑不同边长的小网格内覆盖的点数时，记F中的点落在小方格数目为N(r)，得到拟合ln(N(r))-ln(r)的拟合直线，其中斜率的绝对值就是断裂的容量维值，同样可以得到断裂的拟合度R²。

当考虑不同边长的小网格内覆盖的点数时，统计箱子F中的点落在第i个小网格中的概率，进而可以求取I(r)，得到拟合ln(I(r))-ln(r)的拟合直线，其中斜率的绝对值就是断裂的面积信息维值，同样可以得到断裂的拟合度R²。

本发明的有益效果是：本发明在断裂带数据化后，通过设置边界参数，起始位置，定义统计箱子(F)的参数，通过编写程序使统计箱子(F)不断的移动，并且计算统计箱子(F)所在位置的分维值，实现断裂分形特征的程序化计算。本发明对于断裂分布特征定量评价、断裂分形特征定量描述等多个方面具有较高的实用价值，并且预测成本低廉、可操作性强，大量减少人力的支出，评价结果对预测油气勘探重点区域、裂缝优势分布区域以及低级序断层的分布特征等多个方面有一定的指导意义。

附图说明

图1为一种断裂分形特征优化计算方法的流程图。

图2为断裂带数据化取点次序示意图。

图3为断裂带三点确定其所在的矩形示意图。

图4为按断裂充填步长b对断裂充填示意图。

图5为统计箱子(F)移动轨迹示意图。

图6为金湖凹陷构造位置图。

图7为金湖凹陷阜二段储层断裂分布图。

图8为断裂数据化-充填后示意图。

图9为金湖凹陷阜二段储层断裂容量维分布图。

图10为研究区断裂容量维分形特征曲线。

图11为金湖凹陷阜二段储层断裂面积信息维分布图。

图12为研究区断裂面积信息维分形特征曲线。

具体实施方式

下面结合附图说明本发明的具体实施方式：

金湖凹陷在构造上位于苏北盆地东台坳陷的西部，是晚白垩世发育起来的断陷盆地。北起建湖隆起，南至天长凸起，西邻张八岭隆起，东为菱塘桥和柳堡两个低凸起，地理上跨越江苏、安徽两省，面积约5000km²，是苏北盆地中最大的一个沉积凹陷(图6)。在北东向建湖隆起、天长凸起的控制下，受北东向柳堡低凸起、菱塘桥低凸起及边界杨村断裂的制约，阜宁组沉积时期形成南断北超、南陡北缓的箕状凹陷，并呈北东向展布。

金湖凹陷结构为单断断超式即简单的半地堑，总体表现为东断西超，边界断层为杨村断层。吴堡事件后，北东向断裂发育，单一箕状凹陷在三垛期分为南北两个箕状断陷，北为北东向三河次凹为中心，东为石港断裂带，西为刘庄西北斜坡的东陡西缓的北东向展布的凹陷结构，南为以龙岗次凹为中心，东、西、北三面为斜坡的南陡北缓的箕状凹陷结构，凹陷内自东向西主要发育西部斜坡带、石港断裂带、汉涧斜坡带、卞闵杨构造带、唐湾构造带、宝应斜坡带等正向构造单元和三河次凹、汉涧次凹、龙岗次凹、汜水次凹等四个负向构造单元。

金湖凹陷经历了仪征、吴堡、真武、三垛和盐城等构造运动，多期次构造运动形成了凹陷内部复杂的断裂系统(图7)。金湖凹陷内部发育了不同走向与不同级别的断层。南部的杨村断层作为凹陷主边界断层控制了箕状凹陷的形态，石港断层与铜城断层作为凹陷内二级断层控制了凹陷内次凹的分布。

金湖凹陷阜二段断裂极其复杂，人工计算费时而且计算精度低，断裂的信息维数值通过人工也是无法准确求取的，因此选取该区域对设计的算法测试。首先按照图2所示，对断裂依次取点，完成对该区断裂数字化，统计各类断裂826条，数据化后的点数为22663个。依据区域小断裂的宽度，设定充填步长为10m，对断裂覆盖，按照公式(1)-(3)其中在断裂带内充填的小格子数目为6756359个，设定研究区的起始坐标值为(0m，0m)，终止边界条件为(120000m，120000m)，设定统计箱子(F)为边长为10000m，设定统计箱子(F)的移动步长为1000m，其中，计算了有断裂分布的统计箱子有5781个，将统计箱子依次划分为边长10000m/3、10000m/4、10000m/5、10000m/6、10000m/7、10000m/8(不同r值)的网格，分别计算了它们的信息维值、容量维值，并且分别计算了每个箱子的断裂的相似度(拟合度)R²，其中93％的统计箱子的R²>0.97，作出断裂的分维等值线图以及对应的分形特征曲线(图9-图12所示)。

通过上述计算例，金湖凹陷阜二段储层断裂分布极为复杂，但通过设计的程序计算，断裂分形特征的符合断裂的自相似性，认为本算法的准确可行性，并且运算效率高，能够有效的评价断裂的分形特征，能够有效预测研究区的裂缝分布、断裂活动、油气运聚与断裂分布的关系，为下一步油气的勘探提供指导。

上面以举例方式对本发明进行了说明，但本发明不限于上述具体实施例，凡基于本发明所做的任何改动或变型均属于本发明要求保护的范围。

Claims (4)

1.一种断裂分形特征优化计算方法，所述的步骤如下：

1)根据地震解释构造图，对断裂带数据化，采取折线式数据化方法，分别在断裂的上升盘对应的断层线与下降盘对应的断层线依次取点，其中数据的起点与终点分别对应断裂的两端点，在断层走向变化大的区域，采用密集取点方法，保证取点次序以及断裂的原始形态，将每条断裂读取的点序号标记为(a₀、a₁、a₂...a_n-1、a_n)；

2)设置工区起始位置、充填步长以及边界条件，根据研究区的范围，设置工区的起始范围、边界范围以及断裂带充填步长b(单位m)，断裂带充填步长b主要与断裂带的宽度有关，可设置为工区小断裂宽度的十分之一，或更小；

3)对断裂带充填，利用每条断裂读取点的序号(a₀、a₁、a₂...a_n-1、a_n)，每条断裂由一系列依次排列的三角形组成，对组成断裂带的三角形充填边长为b的小网格；

4)统计箱子(F)参数设置，所述的统计箱子(F)是指在断裂分形中所要计算-评价维数的网格单元，设置统计箱子(F)的边长L(正方形)，统计箱子(F)每次的移动步长为⊿l，通过编程使统计箱子工区逐次移动⊿l；

5)求取断裂的面积信息维与容量维，筛选落入断裂带的所有小格子组成的点集T_sum，通过编写程序，统计落入每个统计箱子(F)中的数目、X坐标、Y坐标，然后将这个统计箱子(F)所在的矩形区域分割成若干个边长为r的正方形后，统计落入每个小网格的数目，通过不断变换r，分别统计落入小方形中小网格数目。

2.根据权利要求1所述的一种断裂分形特征优化计算方法，其特征在于：

所述的对断裂带内充填入边长为b的小网格，原理步骤如下：

对断裂数据化后，每条断裂由一系列依次排列的三角形组成，对于任意的三角形，设其由三点a_m、a_m+1、a_m+2组成，其中m的变化范围0-(n-2)，点a_m、a_m+1、a_m+2对应的X轴、Y轴坐标分别为(x_m，y_m)、(x_m+1，y_m+1)、(x_m+2，y_m+2)，可以依据三个点坐标数值确定充填的范围-矩形ABCD，根据断裂带充填步长b，将矩形ABCD划分为一系列边长为b的小格子(图4)，依据断裂带充填步长b以及矩形ABCD的边界坐标值，可以表示矩形ABCD内任意一个小格子中心所在的坐标数值s(x_p，y_q)，通过公式(1)-(3)，可以判断任意的点s(x_p，y_q)是否落入断裂带(点a_m、a_m+1、a_m+2)组成的三角形内：

参数aa₁、bb₁、cc₁表示为：

$\left\{\begin{array}{ccc}{\mathrm{aa}}_1=(x_{m+2}-x_m)& (y_{m+1}-y_m)-(y_{m+2}-y_m)& (x_{m+1}-x_m)\\ {\mathrm{bb}}_1=(x_{m+1}-x_m)& (x_{m+2}-y_m)-(y_{m+1}-y_m)& (x_{m+2}-x_m)\\ {\mathrm{cc}}_1=(x_m-x_{m+2})& (x_{m+2}-y_{m+1})-(y_m-y_{m+2})& (x_{m+2}-x_{m+1})\end{array}\right.---\left(1\right)$

参数aa₂、bb₂、cc₂表示为：

$\left\{\begin{array}{ccc}{\mathrm{aa}}_1=(x_p-x_m)& (y_{m+1}-y_m)-(y_q-y_m)& (x_{m+1}-x_m)\\ {\mathrm{bb}}_1=(x_p-x_m)& (x_{m+2}-y_m)-(y_q-y_m)& (x_{m+2}-x_m)\\ {\mathrm{cc}}_1=(x_p-x_{m+2})& (x_{m+2}-y_{m+1})-(y_q-y_{m+2})& (x_{m+2}-x_{m+1})\end{array}\right.---\left(2\right)$

参数aa、bb、cc表示为：

$\left\{\begin{array}{c}\mathrm{aa}={\mathrm{aa}}_1\·{\mathrm{aa}}_2\\ \mathrm{bb}={\mathrm{bb}}_1\·{\mathrm{bb}}_2\\ \mathrm{cc}={\mathrm{cc}}_1\·{\mathrm{cc}}_2\end{array}\right.---\left(3\right)$

当aa≥0且bb≥0且cc>0时，认为小格子s(x_p，y_q)落在由点a_m、a_m+1、a_m+2组成的三角形内，依据公式(1)-(3)，在对断裂网格化的基础上，可以筛选出所有落入断裂带的小格子s，在此基础上，实现对断裂带的充填；落在断裂带内小格子s的总数与每个小格子的面积乘积就便可近似为断裂的面积。

3.根据权利要求1所述的一种断裂分形特征优化计算方法，其特征在于：

所述的断裂面积信息维与容量维的求取，其步骤如下：

通过编写程序使统计箱子(F)首先沿X轴移动，每次沿X轴正方向移动步长⊿l，直至到达右边界，完成一次X轴方向移动；之后返回统计箱子(F)起始位置，沿Y轴正方向移动步长⊿l，开始第二次X轴方向移动，直至到达右边界，完成第二次X轴方向移动；之后返回统计箱子(F)起始位置，沿Y轴正方向移动2×⊿l；......直到统计箱子(F)沿Y轴达到上边界时，进行最后一次沿X轴方向移动，到达边界的东北边界，完成在工区范围内的移动。

4.根据权利要求1所述的一种断裂分形特征优化计算方法，其特征在于：

所述的求取断裂的面积信息维与容量维，其步骤如下：

容量维求取方法：当不考虑不同边长的小网格内覆盖的点数时，记录统计箱子(F)中不同边长(r)的小网格的数目为N(r)，变换小网格的边长，得到多个N(r)、r的数值，拟合ln(N(r))-ln(r)的拟合直线，其中斜率的绝对值就是断裂的容量维值，依据拟合的直线，同样可以得到断裂的拟合度R²；

面积信息维求取方法：当考虑不同边长的小网格内覆盖的点数时，统计箱子F中的点落在第i个小网格中的概率，进而可以求取I(r)，变换小网格的边长，得到多个I(r)、r的数值，拟合ln(I(r))-ln(r)的拟合直线，其中斜率的绝对值就是断裂的面积信息维值，依据拟合的直线，同样可以得到断裂的拟合度R²。