CN107578471B - 一种自支撑裂缝初始形态构建方法 - Google Patents

Abstract

本发明公开了一种自支撑裂缝初始形态构建方法，包括：（A）获取目标储层岩石，加工为长方体岩板，将其沿垂直于高度方向一分为二，形成A、B两个粗糙裂缝面；（B）采集A、B两个裂缝面的表面形貌数据，进行降噪和插值处理，得到“初始高度矩阵”A ₀ 、B ₀ ；（C）对裂缝面进行净值化处理，得到“净高度矩阵”A ₁ 、B ₁ ；（D）对裂缝面进行反向处理，得到“反向高度矩阵”A ₂ 、B ₂ ；（E）对裂缝面进行错位处理，得到“错位高度矩阵”A ₃ 、B ₃ ；（F）对裂缝面进行触碰处理，得到“接触高度矩阵”A ₄ 、B ₄ ；（G）利用A ₄ 、B ₄ 作图，得到自支撑裂缝的初始形态。本发明能够构建目标储层自支撑裂缝的初始形态，为清水压裂后产量预测和增产潜力评价提供有效指导。

Description

一种自支撑裂缝初始形态构建方法

技术领域

本发明涉及石油天然气勘探开发领域，尤其是水力压裂领域中一种自支撑裂缝初始形态构建方法。

背景技术

水力压裂技术是低渗透油气藏增产改造的重要措施。水力压裂是利用地面高压泵组，以超过地层吸收能力的排量将压裂液泵入地层来产生裂缝，然后继续注入带有支撑剂(砂粒)的压裂液，使裂缝继续延伸并在其中充填支撑剂，当压裂液返排后，在地层压力作用下，支撑剂在裂缝中起到支撑裂缝的作用，阻止裂缝闭合，从而在地层中形成具有一定长度、允许流体流动的填砂裂缝。

清水压裂是水力压裂的一种形式，被广泛应用于致密砂岩、页岩等油气藏的增产改造中。它具体是指在压裂过程中不加入支撑剂(砂粒)，仅通过将低粘度的压裂液泵入地层来产生裂缝。地下岩石性质差异较大，压裂形成的裂缝表面一般凹凸不平，同时还会在剪切作用下发生错位，因此即使不加入支撑剂，裂缝表面的凸点之间也可以相互支撑形成自支撑裂缝，使裂缝在地层压力作用下仍能保持一定开启程度和流动通道，从而达到改善油气流动条件和油气井增产的目的。由于自支撑裂缝内未充填支撑剂，裂缝形态将直接决定着内部流动通道分布，进而影响油气流量大小。因此，准确构建自支撑裂缝形态对于清水压裂后产量预测和油气藏增产潜力评价具有十分重要的意义。

裂缝初始形态是指裂缝的两个粗糙壁面刚接触时所对应的裂缝形态。自支撑裂缝仅依靠表面凹凸不平的凸起相互支撑以阻止裂缝闭合，但在高地层压力作用下，大量裂缝表面凸起将会受压变形，甚至被压碎，这个过程中裂缝形态变化十分复杂，构建难度也较大。本发明所建立的裂缝形态构建方法则主要是针对自支撑裂缝的初始形态，即自支撑裂缝未发生变形和破碎情况。

噪点在不同的情况下具有不同的含义，在数据采集过程中的噪点主要是指由于仪器本身缺陷或外部干扰而导致结果中出现的一些明显不正常的数据点。

天然裂缝是相对于人工裂缝而言的，指的是地层中由于地壳运动或其他自然因素而天然存在的一类裂缝。天然裂缝一般处于闭合状态，人工裂缝在向前延伸的过程中与天然裂缝相交后，压裂液将侵入天然裂缝，从而促使天然裂缝张开，甚至发生错位。

目前针对自支撑裂缝初始形态构建，国内外学者已提出了多种方法，主要包括随机裂缝形态方法、规则裂缝形态方法以及分形方法，但这些方法普遍存在计算量大、计算过程复杂、构建精度差等缺陷，不具备针对性和实用性，难以进行推广。

发明内容

本发明的目的在于提供一种自支撑裂缝初始形态构建方法，该方法原理可靠，操作简单，能够准确构建出自支撑裂缝的初始形态，进而为清水压裂后产量预测和增产潜力评价提供有效指导。

为达到以上技术目的，本发明提供以下技术方案。

考虑到地下真实自支撑裂缝的初始形态难以直接获取，本发明先获取地层岩石，再运用人工手段将岩石剖开，形成两个粗糙裂缝面并采集得到裂缝表面形貌数据，然后对表面形貌数据进行系列变换以模拟地下真实自支撑裂缝，最后利用变换后的表面形貌数据作图，从而间接构建出自支撑裂缝的初始形态。

一种自支撑裂缝初始形态构建方法，依次包括以下步骤：

(A)获取目标储层岩石，将其加工为长方体岩板，并通过劈裂、剪切或者其它方式将其沿垂直于高度方向一分为二，形成A、B两个粗糙裂缝面；

(B)采集A、B两个裂缝面的表面形貌数据，并对数据进行降噪和插值处理，得到“初始高度矩阵”分别为A₀、B₀；

(C)基于步骤B，对裂缝面进行净值化处理，由A₀、B₀得到“净高度矩阵”A₁、B₁；

(D)基于步骤C，对裂缝面进行反向处理，由A₁、B₁得到“反向高度矩阵”A₂、B₂；

(E)基于步骤D，对裂缝面进行错位处理，由A₂、B₂得到“错位高度矩阵”A₃、B₃；

(F)基于步骤E，对裂缝面进行触碰处理，由A₃、B₃得到“接触高度矩阵”A₄、B₄；

(G)基于步骤F，利用A₄、B₄作图，即可得到自支撑裂缝的初始形态。

所述步骤(B)中，采集的裂缝表面形貌数据主要为裂缝面各点的位置，及其距离岩板底面的垂直距离(即高度)。数据采集可以通过表面轮廓仪实现，一般表面轮廓仪采集得到的数据主要包括裂缝面各点在岩板底面所在平面上的横坐标、纵坐标以及对应高度值。

由于表面轮廓仪本身缺陷或人为干扰等原因，采集得到的表面形貌数据往往存在部分噪点，因此需要进行降噪处理，去除该部分数据点；同时，为满足数据完整性和后续矩阵运算要求，还需要进一步补全噪点区域数据点。常用数据补全方法有很多，本发明优选克里金插值法作为噪点区域的数据补全方法。下面对克里金插值法的操作步骤进行详细介绍(张靖.基于克里金算法的点云数据插值研究[D].长安大学,2014)：

假设采集得到的数据点共有N个，其中噪点有M个，则正常数据点共有(N-M)个。克里金插值主要是根据(N-M)个正常点的坐标、对应高度值和M个噪点的坐标估算M个噪点的高度值。

1)计算所有正常数据点之间的距离和半方差；对于任意两点P_i(i＝1,2,···,N-M)、P_j(j＝1,2,···,N-M且j≠i)，其坐标分别为(x_i,y_i)、(x_j,y_j)，高度值分别为Z(x_i,y_i)、Z(x_j,y_j)，则P_i与P_j之间的距离d_ij和半方差r_ij可通过下式计算：

2)以d_ij为横坐标，r_ij为纵坐标，将d_ij，r_ij绘制成散点图，并进行数据拟合，得到函数关系r＝f(d)；

3)根据函数关系r＝f(d)和所有正常数据点间的距离d_ij重新计算其对应的半方差

4)利用噪点坐标(x_k,y_k)(k＝N-M+1,N-M+2,···,N)计算其与所有正常数据点之间的距离d_ki，并根据函数关系r＝f(d)计算对应的半方差

5)根据

求解以下线性方程组，得到最优权重系数λ_ki；

(其中n＝N-M)

6)根据最优权重系数λ_ki和所有正常数据点的高度值Z(x_i,y_i)，按下式加权得到所有噪点的高度值：

在进行降噪和插值处理后，为了方便后续数据处理，还需将采集得到的数据点按顺序放置于“初始高度矩阵”中，该矩阵的行数和列数取决于表面轮廓仪的数据采集精度。以A裂缝面为例，其对应的“初始高度矩阵”为A₀，则有如下关系存在：

A₀＝(a_ij)_m×n(i＝1,2,…,m；j＝1,2,…,n) (4)

a_ij＝Z(iΔx-Δx,jΔy-Δy) (5)

m＝int(w/Δx)+1 (6)

n＝int(l/Δy)+1 (7)

式中a_ij——矩阵A₀中第i行、第j列元素；

m——矩阵A₀的行数；

n——矩阵A₀的列数；

Z(iΔx-Δx,jΔy-Δy)——横坐标iΔx-Δx、纵坐标jΔy-Δy所对应的裂缝面高度，mm；

Δx——岩板底面所在平面上x方向(对应于岩板宽度方向)的数据采集精度，mm；

Δy——岩板底面所在平面上y方向(对应于岩板长度方向)的数据采集精度，mm；

w——岩板宽度(裂缝宽度)，mm；

l——岩板长度(裂缝长度)，mm；

int(x)——为取整函数，即只截取x中的整数部分，例如：int(5.34)＝5。受数据采集精度影响，w/Δx和l/Δy可能非整数，而数据采集一般是等间距采集，因此该处运用取整函数求取矩阵A₀的行数和列数。

综上所述，可由A裂缝面的表面形貌数据得到其对应的“初始高度矩阵”A₀，B裂缝面的处理方式与A裂缝面相同，即可由B裂缝面的表面形貌数据得到其对应的“初始高度矩阵”B₀。后续对裂缝表面形貌数据的处理则主要是通过矩阵运算的方式进行。

所述步骤(C)中，净值化处理是指用裂缝面上所有点的高度值减去裂缝面最低点的高度值，以此来消除岩板厚度对数据处理结果的影响。净值化处理后将得到A、B两个裂缝面的“净高度矩阵”A₁、B₁，该过程所对应的矩阵运算如下：

A₁＝A₀-[min(A₀)]_m×n (8)

B₁＝B₀-[min(B₀)]_m×n (9)

式中A₁——A裂缝面对应的“净高度矩阵”；

B₁——B裂缝面对应的“净高度矩阵”；

A₀——A裂缝面对应的“初始高度矩阵”；

B₀——B裂缝面对应的“初始高度矩阵”；

min(A₀)——矩阵A₀中最小的元素值；

min(B₀)——矩阵B₀中最小的元素值。

所述步骤(D)中，由于采集裂缝表面形貌数据时两个裂缝面的放置方向是相同的，而在实际裂缝中两个裂缝面却处于相向位置，因此在净值化处理之后还需对其中一个裂缝面进行反向处理(本发明以A裂缝面为例处理)。同时，为了保证反向处理后两个裂缝面不发生重叠，还需要将反向后的裂缝面沿高度方向整体平移一定距离。反向处理后将得到A、B两个裂缝面的“反向高度矩阵”A₂、B₂，该过程所对应的矩阵运算如下：

A₂＝[max(A₁)]_m×n-A₁+[max(B₁)]_m×n (10)

B₂＝B₁ (11)

式中A₂——A裂缝面对应的“反向高度矩阵”；

B₂——B裂缝面对应的“反向高度矩阵”；

A₁——A裂缝面对应的“净高度矩阵”；

B₁——B裂缝面对应的“净高度矩阵”；

max(A₁)——矩阵A₁中最大的元素值；

max(B₁)——矩阵B₁中最大的元素值。

所述步骤(E)中，在地层条件下裂缝的两个粗糙面会发生错位，从而依靠裂缝表面凸起相互支撑形成自支撑裂缝，因此在反向处理后需要对两个裂缝面进一步进行错位处理。该过程的主要思想是删除A裂缝面长度方向末端和B裂缝面长度方向起始端的数据点，具体删除范围取决于裂缝面的错位量。错位处理后将得到A、B两个裂缝面的“错位高度矩阵”A₃、B₃，同时由于删除了部分数据点，矩阵A₃、B₃的行数仍与A₁、B₁、A₂、B₂保持一致，但列数将会减少。该过程所对应的矩阵运算如下：

A₃＝(b_ij)_m×n′(i＝1,2,…,m；j＝1,2,…,n′) (12)

B₃＝(c_ij)_m×n′(i＝1,2,…,m；j＝u+1,u+2,…,n) (13)

n′＝n-int(δ/Δy)＝int(l/Δy)-int(δ/Δy)+1 (14)

u＝int(δ/Δy) (15)

式中A₃——A裂缝面对应的“错位高度矩阵”；

B₃——B裂缝面对应的“错位高度矩阵”；

b_ij——A裂缝面对应的“反向高度矩阵”A₂中第i行、第j列元素；

c_ij——B裂缝面对应的“反向高度矩阵”B₂中第i行、第j列元素；

m——矩阵A₃、B₃的行数，同矩阵A₁、B₁、A₂、B₂；

n′——矩阵A₃、B₃的列数；

n——矩阵A₂、B₂的列数，同矩阵A₁、B₁；

int(x)——取整函数，同上；

δ——裂缝错位量，mm；

Δy——岩板底面所在平面上y方向(对应于岩板长度方向)的数据采集精度，mm；

u——矩阵A₂、B₂中需删除的数据列数。

上述式子中，式(14)、(15)中的裂缝错位量δ是一个与地层所处应力环境和地层岩石力学性质有关的未知量，可通过下面的公式求取：

式中δ——裂缝错位量；

E——地层岩石的杨氏模量，MPa；

v——地层岩石的泊松比，无量纲；

L——人工裂缝半长，m；

Δσ——地层水平主应力差，Δσ＝σ₃-σ₁，MPa；

σ₃——地层最大水平主应力，MPa；

σ₁——地层最小水平主应力，MPa；

θ——人工裂缝与地层最大水平主应力方向之间的夹角，度。

该公式的推导过程如下：

1)假设地层中存在一条长度为2L₀的闭合天然裂缝，其受地层最大水平主应力σ₃和最小水平主应力σ₁作用，且与最大水平主应力方向的夹角为θ₀(0＜θ₀＜π/2)。根据线弹性理论，作用在裂缝面上的法向应力和剪切应力分别为(周健,陈勉,金衍,等.压裂中天然裂缝剪切破坏机制研究[J].岩石力学与工程学报,2008,27:2637—2641)：

式中σ——最大、最小水平主应力在天然裂缝面上产生的法向应力，MPa；

τ——最大、最小水平主应力在天然裂缝面上产生的剪切应力，MPa；

σ₃——地层最大水平主应力，MPa；

σ₁——地层最小水平主应力，MPa；

θ₀——天然裂缝与地层最大水平主应力方向的夹角，度。

2)地层中一般发育有大量的天然裂缝，在压裂过程中，人工裂缝与天然裂缝相交前，作用于天然裂缝面上的法向应力为压应力，裂缝处于闭合状态。当人工裂缝与天然裂缝相交后，随着压裂液侵入天然裂缝，天然裂缝面上的法向应力减小，同时较大的剪切应力将促使裂缝面发生剪切错位。天然裂缝面上的有效剪切应力可由裂缝面上的法向应力σ和剪切应力τ表示(Warpinski N R,Teufel L W.Influence of geologic discontinuities onhydraulic fracture propagation.Journal of Petroleum Technology,SPE13224,1987:209—220)：

τ_e＝|τ|-τ₀-K_f(σ-P) (19)

式中τ_e——天然裂缝面上有效剪切应力，MPa；

τ——最大、最小水平主应力在天然裂缝面上产生的剪切应力，MPa；

σ——最大、最小水平主应力在天然裂缝面上产生的法向应力，MPa；

τ₀——天然裂缝面上的粘聚力，为天然裂缝的固有性质，MPa；

K_f——天然裂缝面上的摩擦系数，为天然裂缝的固有性质，MPa；

P——天然裂缝中流体压力，MPa。

3)天然裂缝面在有效剪切应力τ_e作用下发生的错位量可通过断裂力学理论获得。根据断裂力学中的Westergaard函数，在无限大介质中裂缝面任意一点的剪切错位量可表示为(李世愚,和泰名,尹祥础,等.岩石断裂力学导论.合肥:中国科学技术大学出版社,2010:50—51)：

式中ω(x)——裂缝长度方向坐标x处的错位量，mm；

k——Kolosov常数，k＝3-4v；

G——地层岩石的剪切模量，

E——地层岩石的杨氏模量，MPa；

v——地层岩石的泊松比，无量纲；

τ_e——裂缝面上的有效剪切应力，MPa；

x——裂缝长度方向上任意点坐标，m；

L₀——天然裂缝半长，m。

由式(20)可知，在裂缝中心位置处(x＝0)错位量可取得最大值δ，具体表示为：

式中δ——裂缝面最大错位量；

E——地层岩石的杨氏模量，MPa；

v——地层岩石的泊松比，无量纲；

L₀——天然裂缝半长，m；

τ_e——裂缝面上的有效剪切应力，MPa。

4)在压裂过程中流体压力会不断发生波动，同时天然裂缝面上的粘聚力以及摩擦阻力均属于不确定因素(不同位置具有不同的值)。为减弱上述因素对裂缝面错位量计算结果的影响，本发明假设天然裂缝处于完全张开状态(P＝σ)，忽略天然裂缝面上的摩擦阻力(K_f＝0)和粘聚力(τ₀＝0)，并将裂缝中心位置处的错位量作为裂缝实际错位量，则通过式(17)、(18)、(19)、(21)可获得裂缝面的错位量表达式如下：

式中δ——裂缝面错位量；

E——地层岩石的杨氏模量，MPa；

v——地层岩石的泊松比，无量纲；

L₀——天然裂缝半长，m；

Δσ——地层水平主应力差，Δσ＝σ₃-σ₁，MPa；

σ₃——地层最大水平主应力，MPa；

σ₁——地层最小水平主应力，MPa；

θ₀——天然裂缝与地层最大水平主应力方向之间的夹角，度。

特别的，虽然上述推导过程是针对地层中的天然裂缝，但由于人工裂缝与天然裂缝处于相同的地应力环境，具有相似的岩石力学性质，同时式(22)并未考虑天然裂缝中摩擦阻力、粘聚力等特有属性，因此可以用式(22)估算人工裂缝的错位量，只需将式(22)中的L₀用人工裂缝半长L代替、θ₀用人工裂缝与地层最大水平主应力方向之间的夹角θ代替即可。

所述步骤(F)中，在初始状态下，自支撑裂缝的两个面是相互接触的，因此在错位处理后还需进一步对裂缝进行触碰处理。触碰处理是指不断缩小两个裂缝面之间的间距，直至裂缝面刚好接触。因此，可先求取错位处理后裂缝A面与B面之间的最小高度差h，然后将A裂缝面向下平移h，即可实现两个裂缝面的触碰。触碰处理后将得到A、B两个裂缝面的“接触高度矩阵”A₄、B₄，该过程所对应的矩阵运算如下：

h＝min(d_ij-e_ij)(i＝1,2,…,m；j＝1,2,…,n′) (23)

A₄＝A₃-[h]_m×n′ (24)

B₄＝B₃ (25)

式中h——错位处理后A裂缝面与B裂缝面之间的最小高度差，mm；

d_ij——A裂缝面对应的“错位高度矩阵”A₃中第i行、第j列元素；

e_ij——B裂缝面对应的“错位高度矩阵”B₃中第i行、第j列元素；

A₄——A裂缝面对应的“接触高度矩阵”；

B₄——B裂缝面对应的“接触高度矩阵”；

m——“错位高度矩阵”A₃、B₃的行数；

n′——“错位高度矩阵”A₃、B₃的列数。

所述步骤(G)中，由于“接触高度矩阵”A₄、B₄中仅保留着A、B两个裂缝面经过一系列处理之后的高度值，并没有具体对应的横、纵坐标，因此需要先通过下面关系式将“接触高度矩阵”A₄、B₄中的高度值转化为空间中裂缝面各个坐标点对应的高度值，再通过作图软件或其它方式进行作图，最终得到自支撑裂缝的初始形态。

Z_A(iΔx-Δx,jΔy-Δy)＝f_ij(i＝1,2,…,m；j＝1,2,…,n′) (26)

Z_B(iΔx-Δx,jΔy-Δy)＝g_ij(i＝1,2,…,m；j＝1,2,…，n′) (27)

式中Z_A(iΔx-Δx,jΔy-Δy)——A裂缝面经一系列处理后，横坐标为iΔx-Δx，纵坐标为jΔy-Δy所对应点的高度，mm；

Z_B(iΔx-Δx,jΔy-Δy)——B裂缝面经一系列处理后，横坐标为iΔx-Δx，纵坐标为jΔy-Δy所对应点的高度，mm；

f_ij——A裂缝面对应的“接触高度矩阵”A₄中第i行、第j列元素；

g_ij——B裂缝面对应的“接触高度矩阵”B₄中第i行、第j列元素；

m——“接触高度矩阵”A₄、B₄的行数；

n′——“接触高度矩阵”A₄、B₄的列数。

与现有技术相比，本发明的有益效果是：提出了一种自支撑裂缝初始形态的构建方法，该方法原理可靠，操作简便，能够准确构建目标储层自支撑裂缝的初始形态，从而为目标储层清水压裂后产量预测和增产潜力评价提供有效指导。

附图说明

图1是对裂缝面进行净值化处理示意图(以x＝18mm截面为例)。

图2是对裂缝面进行反向处理示意图(以x＝18mm截面为例)。

图3是对裂缝面进行错位处理示意图(以x＝18mm截面为例)。

图4是对裂缝面进行触碰处理示意图(以x＝18mm截面为例)。

图5是川南“志留系龙马溪组”储层自支撑裂缝初始形态。

具体实施方式

下面结合附图和实施例来进一步详细说明本发明。

川南页岩气是西南探区页岩气勘探开发先导实验区，其中“志留系龙马溪组”是主要的页岩气储集层，也是增产改造的目标层位。该区块前期已进行大量清水压裂施工，并取得了一定成果，但对于压裂后在地层中所形成的裂缝形态仍然缺乏有效构建方法，因此开发人员往往难以准确预测压裂后的产量，进而评价压裂后的增产效果。

以川南页岩气“志留系龙马溪组”储层为例，采用本发明所提供的方法进行自支撑裂缝初始形态构建，主要步骤如下：

(A)获取川南页岩气“志留系龙马溪组”储层岩石，将其加工为长178mm、宽36mm、高50mm的长方体岩板，并通过劈裂方式将其沿垂直于高度方向一分为二，形成A、B两个粗糙裂缝面；

(B)通过表面轮廓仪采集A、B两个裂缝面的表面形貌数据，表面轮廓仪在x(岩板宽度方向)和y方向(岩板长度方向)的采集精度均为0.02mm；然后对数据进行降噪和插值处理，得到“初始高度矩阵”分别为A₀、B₀。由于采集得到的数据点较多，相对应的“初始高度矩阵”规模也较大(1800行，8900列)，受篇幅所限，此处未详细列出矩阵A₀、B₀中各元素的具体数值，下同；

(C)对裂缝面进行净值化处理(见图1)，由A₀、B₀得到“净高度矩阵”A₁、B₁；

(D)对裂缝面进行反向处理(见图2)，由A₁、B₁得到“反向高度矩阵”A₂、B₂；

(E)对裂缝面进行错位处理(见图3)，由A₂、B₂得到“错位高度矩阵”A₃、B₃；该过程中错位量的计算过程如下：

川南“志留系龙马溪”组页岩最大水平主应力σ₃为89.37MPa，最小水平主应力σ₁为83.6MPa，水平主应力差Δσ为5.77MPa，平均杨氏模量为23460MPa，平均泊松比为0.237。根据前期微地震监测经验，人工裂缝长度为20m左右，与最大水平主应力方向夹角为45°左右，则利用式(16)可计算得到错位量δ＝4.6mm；

(F)对裂缝面进行触碰处理(见图4)，由A₃、B₃得到“接触高度矩阵”A₄、B₄；

(G)基于“接触高度矩阵”A₄、B₄，利用式(26)、(27)得到经过一系列处理后A、B裂缝面上各点所处的空间位置，并通过作图软件进行作图，最终得到“志留系龙马溪组”储层中自支撑裂缝初始形态(见图5)。

Claims (2)

1.一种自支撑裂缝初始形态构建方法，依次包括以下步骤：

(A)获取目标储层岩石，加工为长方体岩板，将其沿垂直于高度方向一分为二，形成A、B两个粗糙裂缝面；

(B)采集A、B两个裂缝面的表面形貌数据，对数据进行降噪和插值处理，得到“初始高度矩阵”分别为A₀、B₀，包括：将采集得到的表面形貌数据进行降噪和插值处理后，对于A裂缝面，得到其对应的“初始高度矩阵”A₀：

A₀＝(a_ij)_m×n，i＝1,2,…,m；j＝1,2,…,n

a_ij＝Z(iΔx-Δx,jΔy-Δy)

m＝int(w/Δx)+1

n＝int(l/Δy)+1

式中a_ij——矩阵A₀中第i行、第j列元素；

m——矩阵A₀的行数；

n——矩阵A₀的列数；

Z(iΔx-Δx,jΔy-Δy)——横坐标iΔx-Δx、纵坐标jΔy-Δy所对应的裂缝面高度，mm；

Δx——岩板底面所在平面上x方向的数据采集精度，mm；

Δy——岩板底面所在平面上y方向的数据采集精度，mm；

w——岩板宽度，即裂缝宽度，mm；

l——岩板长度，即裂缝长度，mm；

int(x)——为取整函数，即只截取x中的整数部分；

对于B裂缝面，得到其对应的“初始高度矩阵”B₀；

(C)对裂缝面进行净值化处理，由A₀、B₀得到“净高度矩阵”A₁、B₁，是指用裂缝面上所有点的高度值减去裂缝面最低点的高度值，以此来消除岩板厚度对数据处理结果的影响，净值化处理后得到A、B两个裂缝面的“净高度矩阵”A₁、B₁：

A₁＝A₀-[min(A₀)]_m×n

B₁＝B₀-[min(B₀)]_m×n

式中A₁——A裂缝面对应的“净高度矩阵”；

B₁——B裂缝面对应的“净高度矩阵”；

min(A₀)——矩阵A₀中最小的元素值；

min(B₀)——矩阵B₀中最小的元素值；

(D)对A裂缝面进行反向处理，由A₁、B₁得到“反向高度矩阵”A₂、B₂：

A₂＝[max(A₁)]_m×n-A₁+[max(B₁)]_m×n

B₂＝B₁

式中A₂——A裂缝面对应的“反向高度矩阵”；

B₂——B裂缝面对应的“反向高度矩阵”；

max(A₁)——矩阵A₁中最大的元素值；

max(B₁)——矩阵B₁中最大的元素值；

(E)删除A裂缝面长度方向末端和B裂缝面长度方向起始端的数据点，由A₂、B₂得到“错位高度矩阵”A₃、B₃：

A₃＝(b_ij)_m×n′，i＝1,2,…,m；j＝1,2,…,n′

B₃＝(c_ij)_m×n′，i＝1,2,…,m；j＝u+1,u+2,…,n

n′＝n-int(δ/Δy)＝int(l/Δy)-int(δ/Δy)+1

u＝int(δ/Δy)

式中A₃——A裂缝面对应的“错位高度矩阵”；

B₃——B裂缝面对应的“错位高度矩阵”；

b_ij——A裂缝面对应的“反向高度矩阵”A₂中第i行、第j列元素；

c_ij——B裂缝面对应的“反向高度矩阵”B₂中第i行、第j列元素；

m——矩阵A₃、B₃的行数，同矩阵A₁、B₁、A₂、B₂；

n′——矩阵A₃、B₃的列数；

n——矩阵A₂、B₂的列数，同矩阵A₁、B₁；

int(x)——取整函数，同上；

δ——裂缝错位量，mm；

u——矩阵A₂、B₂中需删除的数据列数；

(F)不断缩小两个裂缝面之间的间距，直至裂缝面刚好接触，先求取错位处理后裂缝A面与B面之间的最小高度差h，然后将A裂缝面向下平移h，实现两个裂缝面的触碰，由A₃、B₃得到“接触高度矩阵”A₄、B₄：

h＝min(d_ij-e_ij)，i＝1,2,…,m；j＝1,2,…,n′

A₄＝A₃-[h]_m×n′

B₄＝B₃

式中h——错位处理后A裂缝面与B裂缝面之间的最小高度差，mm；

d_ij——A裂缝面对应的“错位高度矩阵”A₃中第i行、第j列元素；

e_ij——B裂缝面对应的“错位高度矩阵”B₃中第i行、第j列元素；

A₄——A裂缝面对应的“接触高度矩阵”；

B₄——B裂缝面对应的“接触高度矩阵”；

m——“错位高度矩阵”A₃、B₃的行数；

n′——“错位高度矩阵”A₃、B₃的列数；

(G)通过下式将“接触高度矩阵”A₄、B₄中的高度值转化为空间中裂缝面各个坐标点对应的高度值，再利用A₄、B₄作图，即可得到自支撑裂缝的初始形态：

Z_A(iΔx-Δx,jΔy-Δy)＝f_ij，i＝1,2,…,m；j＝1,2,…,n′

Z_B(iΔx-Δx,jΔy-Δy)＝g_ij，i＝1,2,…,m；j＝1,2,…,n′

式中Z_A(iΔx-Δx,jΔy-Δy)——A裂缝面经一系列处理后，横坐标为iΔx-Δx，纵坐标为jΔy-Δy所对应点的高度，mm；

Z_B(iΔx-Δx,jΔy-Δy)——B裂缝面经一系列处理后，横坐标为iΔx-Δx，纵坐标为jΔy-Δy所对应点的高度，mm；

f_ij——A裂缝面对应的“接触高度矩阵”A₄中第i行、第j列元素；

g_ij——B裂缝面对应的“接触高度矩阵”B₄中第i行、第j列元素；

m——“接触高度矩阵”A₄、B₄的行数；

n′——“接触高度矩阵”A₄、B₄的列数。

2.如权利要求1所述的一种自支撑裂缝初始形态构建方法，其特征在于，所述裂缝错位量δ通过下式求取：

式中δ——裂缝错位量；

E——地层岩石的杨氏模量，MPa；

v——地层岩石的泊松比，无量纲；

L——人工裂缝半长，m；

Δσ——地层水平主应力差，Δσ＝σ₃-σ₁，MPa；

σ₃——地层最大水平主应力，MPa；

σ₁——地层最小水平主应力，MPa；

θ——人工裂缝与地层最大水平主应力方向之间的夹角，度。

Images