水力压裂改造是实现页岩气井增产的主要手段, 根据水力压裂微地震数据估计矩张量对于研究地震破裂过程的细节、理解非常规储层的地质力学特征具有重要意义^[1-2]。微地震震源机制是指微地震发生的物理力学过程, 通过微地震震源机制, 可以深入分析发震的内/外在诱因以及地下岩体的破(断)裂机理, 研究水力压裂诱发和增强的裂纹方向及裂纹网络的复杂程度等。由于页岩等致密储层岩石具有各向异性^[3], 因此需要考虑各向异性对于微地震震源的影响。微地震矩张量中包含非双力偶分量, 产生的原因多种多样, 各向异性会明显地影响矩张量中的非双力偶分量^[4-5]。VAVRYČUK^[6]指出, 当震源区介质是各向同性时, 如果孔隙压力较小, 可能仅仅产生剪切破裂, 震源机制主要是双力偶类型, 而当裂纹中的孔隙压力较大时, 会发生张破裂, 震源机制中会包含非双力偶分量。RUTLEDGE等^[7]指出, 震源区的各向异性会使得地震源特征变得复杂。VAVRYČUK^[8]研究了各向异性介质下剪切位错源的地震矩张量, 并从数学上确认了各向异性源本身是造成非双力偶机制的一个重要原因。蔡晓刚等^[9]研究了各向异性ATI介质中的剪切位错源地震矩张量。水力压裂区作为应力较为集中的区域, 包含天然裂缝以及水力压裂裂缝, 在整体上呈现各向异性特征^[10], 需要就水力压裂区的各向异性参数对震源机制和矩张量的影响进行深入研究和分析。本文给出了HTI介质中地震矩张量的各向异性解析表达式, 研究了HTI各向异性对双力偶分量(DC)、各向同性分量(ISO)和补偿线性偶极子分量(CLVD)的影响。

1 各向异性HTI介质中的剪张源模型 1.1 含流体HTI介质参数分析

KRIEF等^[11]认为, 背景各向同性岩石的体积模量K与剪切模量μ和固体基质的体积模量K_m与剪切模量μ_m满足以下关系:

$ K = {K_{\rm{m}}}{\left( {1 - \varphi } \right)^{3/\left( {1 - \varphi } \right)}} $

(1)

$ \mu = {\mu _{\rm{m}}}{\left( {1 - \varphi } \right)^{3/\left( {1 - \varphi } \right)}} $

(2)

式中:φ是孔隙度。公式(1)和公式(2)是经验速度-孔隙度公式, 能够模拟不同孔隙度的岩石弹性特征, 即背景各向同性介质的弹性模量参数。

SCHOENBERG等^[12]将裂缝当成一个特殊界面来处理, 应力通过该界面时是连续的, 而位移存在不连续现象。对于含垂直裂纹的HTI介质(图 1), 可以利用该方法获得干燥岩石的刚度矩阵表达式:

$ {C^{{\rm{dry}}}} = \left( {\begin{array}{*{20}{c}} {Q\left( {1 - {\Delta _{\rm{N}}}} \right)}&{\lambda \left( {1 - {\Delta _{\rm{N}}}} \right)}&{\lambda \left( {1 - {\Delta _{\rm{N}}}} \right)}&0&0&0\\ {\lambda \left( {1 - {\Delta _{\rm{N}}}} \right)}&{Q\left( {1 - {r^2}{\Delta _{\rm{N}}}} \right)}&{\lambda \left( {1 - r{\Delta _{\rm{N}}}} \right)}&0&0&0\\ {\lambda \left( {1 - {\Delta _{\rm{N}}}} \right)}&{\lambda \left( {1 - r{\Delta _{\rm{N}}}} \right)}&{Q\left( {1 - {r^2}{\Delta _{\rm{N}}}} \right)}&0&0&0\\ 0&0&0&\mu &0&0\\ 0&0&0&0&{\mu \left( {1 - {\Delta _{\rm{T}}}} \right)}&0\\ 0&0&0&0&0&{\mu \left( {1 - {\Delta _{\rm{T}}}} \right)} \end{array}} \right) $

(3)

式中:Q是纵波模量, Q=λ+2μ, λ是背景各向同性介质的拉梅参数, 满足λ=K－2μ/3;Δ_N和Δ_T分别是法向和切向裂隙弱度, 满足Δ_N=Z_NQ/(1+Z_NQ), Δ_T=Z_Tμ/(1+Z_Tμ), Z_N是附加的裂隙法向韧度, Z_T是附加的裂隙切向韧度; r=λ/(λ+2μ)。

GASSMANN将各向异性岩石的饱和刚度矩阵表述为干燥岩石刚度矩阵和一个附加的流体影响项, 因此利用各向异性GASSMANN方程可以获得饱和各向异性岩石的刚度矩阵各分量^[13-14]:

$ \begin{array}{l} C_{11}^{{\rm{sat}}} = Q\left( {1 - {\Delta _{\rm{N}}}} \right) + \\ \frac{{{{\left[ {{K_{\rm{m}}} - K\left( {1 - {\Delta _{\rm{N}}}} \right)} \right]}^2}}}{{\left( {{K_{\rm{m}}}/{K_{{\rm{fl}}}}} \right)\varphi \left( {{K_{\rm{m}}} - {K_{{\rm{fl}}}}} \right) + {K_{\rm{m}}} - K + {K^2}{\Delta _{\rm{N}}}/Q}} \end{array} $

(4)

$ \begin{array}{l} C_{33}^{{\rm{sat}}} = Q\left( {1 - {r^2}{\Delta _{\rm{N}}}} \right) + \\ \frac{{{{\left[ {{K_{\rm{m}}} - K\left( {1 - r{\Delta _{\rm{N}}}} \right)} \right]}^2}}}{{\left( {{K_{\rm{m}}}/{K_{{\rm{fl}}}}} \right)\varphi \left( {{K_{\rm{m}}} - {K_{{\rm{fl}}}}} \right) + {K_{\rm{m}}} - K + {K^2}{\Delta _{\rm{N}}}/Q}} \end{array} $

(5)

$ \begin{array}{l} C_{13}^{{\rm{sat}}} = \lambda \left( {1 - {\Delta _{\rm{N}}}} \right) + \\ \frac{{\left[ {{K_{\rm{m}}} - K\left( {1 - {\Delta _{\rm{N}}}} \right)} \right]\left[ {{K_{\rm{m}}} - K\left( {1 - r{\Delta _{\rm{N}}}} \right)} \right]}}{{\left( {{K_{\rm{m}}}/{K_{{\rm{fl}}}}} \right)\varphi \left( {{K_{\rm{m}}} - {K_{{\rm{fl}}}}} \right) + {K_{\rm{m}}} - K + {K^2}{\Delta _{\rm{N}}}/Q}} \end{array} $

(6)

$ C_{44}^{{\rm{sat}}} = \mu $

(7)

$ C_{55}^{{\rm{sat}}} = \mu \left( {1 - {\Delta _{\rm{T}}}} \right) $

(8)

式中:K_fl是孔隙流体的体积模量。

本文采用的模型基质体积模量满足K_m=37GPa, 基质剪切模量μ_m=44GPa, 基质密度为2.65g/cm³。首先采用KRIEF等^[11]给出的方法加入孔隙形成各向同性介质, 然后将裂隙参数加入各向同性背景中形成干燥HTI介质。选择水作为孔隙流体, 流体密度ρ_fl=1.02g/cm³, 流体体积模量K_fl=2.7GPa, 在孔隙和裂缝中加入流体形成饱和HTI介质。图 2给出了孔隙度φ=0.2, Δ_T=0.2时干燥、饱水HTI介质的参数变化。由图 2可知, 不同的流体饱和状态会影响纵波速度分布、横波分裂程度和快横波偏振方向, 根据横波分裂的观测结果可以反演获得地下裂隙的状态。

1.2 HTI介质中剪张源的源张量与矩张量

假定裂纹面的法向方向为n=(n₁, n₂, n₃)^T, 位错方向为v=(v₁, v₂, v₃)^T, 图 3为剪张源模型的示意图, 采用走向角φ、倾向角δ、滑动角θ和张裂角α表示的裂纹面参数满足^[15]:

$ {n_1} = - \sin \delta \sin \varphi ,{n_2} = \sin \delta \cos \varphi ,{n_3} = - \cos \delta $

(9)

$ \begin{array}{l} {v_1} = \left( {\cos \theta \cos \varphi + \cos \delta \sin \theta \sin \varphi } \right)\cos \alpha - \\ \;\;\;\;\;\sin \delta \sin \varphi \sin \alpha \\ {v_2} = \left( {\cos \theta \sin \varphi - \cos \delta \sin \theta \cos \varphi } \right)\cos \alpha + \\ \;\;\;\;\;\sin \delta \cos \varphi \sin \alpha \\ {v_3} = - \sin \theta \sin \delta \cos \alpha - \cos \delta \sin \alpha \end{array} $

(10)

走向角φ、倾向角δ、滑动角θ和张裂角α通常是空间位置和时间的函数, 为了简化可以忽略其时空变化, 将剪张错动近似为点源。走向角、倾向角和滑动角描述了断层面的法向和沿着断层面位错矢量的切向, 而张裂角描述了偏离断层面的位错矢量的偏离程度^[16-17]。张裂角α∈(－90°, 90°), 当α > 0时代表剪张源, 对于纯张裂源, α=90°; α < 0代表压裂源, 对于纯压裂源, α=-90°; 对于纯剪切源而言, α=0。

采用裂纹面的法向和位错方向可以获得剪张源的源张量D:

$ \mathit{\boldsymbol{D}} = \frac{{uS}}{2}\left[ {\begin{array}{*{20}{c}} {2{n_1}{v_1}}&{{n_1}{v_2} + {n_2}{v_1}}&{{n_1}{v_3} + {n_3}{v_1}}\\ {{n_1}{v_2} + {n_2}{v_1}}&{2{n_2}{v_2}}&{{n_2}{v_3} + {n_3}{v_2}}\\ {{n_1}{v_3} + {n_3}{v_1}}&{{n_2}{v_3} + {n_3}{v_2}}&{2{n_3}{v_3}} \end{array}} \right] $

(11)

式中:u是位错矢量的大小, S为断层面积, uS为源强度。

各向异性介质微地震矩张量可以由源张量和震源区介质参数获得:

$ {\mathit{\boldsymbol{M}}_{pq}} = {\mathit{\boldsymbol{C}}_{pqkl}} + {\mathit{\boldsymbol{D}}_{kl}} $

(12)

式中:C_pqkl(p, q, k, l=1, 2, 3)为震源区的各向异性刚度参数, 四阶弹性张量C_pqkl可以与二阶弹性张量C_ij(i, j=1, 2, 3, 4, 5, 6)相互转化。

2 各向异性HTI介质中的震源特征分析

采用Hudson图可以对微地震震源机制进行分类描述^[18]。Hudson图中横坐标表示剪切组分T值, 纵坐标表示张性组分k值, T, k的取值范围均为-1~1。此外还可以采用钻石图来表示震源机制中的各个组分^[4]。

2.1 矩张量分解

釆用特征值分解法可以将微地震震源矩张量分解为双力偶部分(DC)、补偿线性矢量极偶成分(CLVD)以及各向同性部分(ISO)^[19], 求取地震矩张量M的特征值和特征方向。设特征值为M₁, M₂和M₃, |M₁|≤|M₂|≤|M₃|, 则各向同性分量M^ISO=(M₁+M₂+M₃)/3。

$ \left[ {\begin{array}{*{20}{c}} {{M_1}}\\ {{M_2}}\\ {{M_3}} \end{array}} \right] = {M^{{\rm{ISO}}}}\left[ {\begin{array}{*{20}{c}} 1\\ 1\\ 1 \end{array}} \right] + \left[ {\begin{array}{*{20}{c}} {M_1^ * }\\ {M_2^ * }\\ {M_3^ * } \end{array}} \right] $

(13)

其中第二项表示偏量部分(DC与CLVD分量之和)M^*, 有多种分解结果, 这里分析DC和CLVD的分解。设|M₁^*|≤|M₂^*|≤|M₃^*|, 则:

$ \left[ {\begin{array}{*{20}{c}} {M_1^ * }\\ {M_2^ * }\\ {M_3^ * } \end{array}} \right] = {M^{{\rm{DC}}}}\left[ {\begin{array}{*{20}{c}} 0\\ { - 1}\\ 1 \end{array}} \right] + \frac{1}{2}{M^{{\rm{CLVD}}}}\left[ {\begin{array}{*{20}{c}} { - 1}\\ { - 1}\\ 2 \end{array}} \right] $

(14)

其中M^DC=M₁^*－M₂^*, M^CLVD=－2M₁^*。采用参数ε描述偏量部分偏离双力偶分量的程度, 则:

$ \varepsilon = \frac{{M_1^ * }}{{M_3^ * }} = \frac{{{M^{{\rm{CLVD}}}}}}{{2\left| {{M^{{\rm{CLVD}}}} + {M^{{\rm{DC}}}}} \right|}} $

(15)

对于纯CLVD源, ε=±0.5;对于纯DC源, ε=0;对于膨胀源, ε是正的, 对于收缩源, ε是负的。

ISO, DC和CLVD分量在地震矩张量中所占的比例满足^[20]:

$ \begin{array}{l} {H_{{\rm{ISO}}}} = \frac{1}{3}\frac{{{M^{{\rm{ISO}}}}}}{{\left| {{M_3}} \right|}} \times 100\% \\ {H_{{\rm{CLVD}}}} = 2\varepsilon \left( {100\% - \left| {{H_{{\rm{ISO}}}}} \right|} \right)\\ {H_{{\rm{DC}}}} = 100\% - \left| {{H_{{\rm{ISO}}}}} \right| - \left| {{H_{{\rm{CLVD}}}}} \right| \end{array} $

(16)

其中H_ISO和H_CLVD可正可负, H_ISO > 0代表体积膨胀, 相反H_ISO < 0代表体积压缩; H_CLVD的正负由ε控制; H_DC为正。

2.2 非DC分量变化特征分析

各向异性介质中的矩张量参数依赖于各向异性的类型和强度以及断层面的方向等参数, 图 4为HTI裂隙介质中由于剪切断裂产生的非DC分量, 其中断层面走向30°, 倾向角为40°, 滑动角为20°, 张裂角从-90°变为90°。图 4a为Hudson图, 其中黑色表示各向同性介质, 红色表示干燥介质, 绿色表示饱和HTI介质; 图 4b为各向同性介质中的矩张量分解结果; 图 4c为干燥HTI介质中的矩张量分解结果; 图 4d为饱水HTI介质中的矩张量分解结果。由图 4可知饱和介质中的ISO分量和DC分量相对干燥介质有所增加。此外, 对于各向同性介质, 当张裂角为0时, 为纯剪切源, 只包含DC分量; 而介质为各向异性时, 由于震源区介质参数的影响, 矩张量中仍然包含非DC分量。研究表明, 较小角度的张裂角也会产生较大的非DC分量。

图 5给出了背景介质孔隙度变化时的剪张源矩张量分解结果, 背景介质孔隙度从0变化到0.3, 其中断层面走向30°, 倾向角为40°, 滑动角为20°, 张裂角为10°。可以看出, 背景介质孔隙度变化会对剪张源矩张量产生影响。各向同性介质中各分量比例不随孔隙度的变化而变化, 而各向异性介质由于背景介质孔隙度的变化会影响各向异性参数, 因而会对震源矩张量产生影响。

图 6给出了裂隙切向弱度参数变化时的剪张源矩张量分解结果, 裂隙切向弱度参数从0变化到0.3, 剪张源参数同图 5。可以看出, 裂隙切向弱度参数变化会对剪张源矩张量产生影响, 震源区各向异性介质由于切向弱度参数的变化会影响各向异性参数, 因而会对震源矩张量产生影响。

3 震源沙滩球与辐射花样分析 3.1 矩张量变化特征分析

图 7给出了不同介质中不同张裂角条件下的震源机制沙滩球, 从上到下依次为各向同性介质、干燥HTI各向异性介质和饱和HTI各向异性介质, 从左到右张裂角分别为0, 10°, 50°, 90°。从图 7可知背景介质会对震源沙滩球中的极性分布产生影响。

3.2 各向异性介质中张裂源的辐射花样

依据HTI介质各向异性地震矩张量解析表达式, 结合Aki和Richards公式, 计算了远场P波辐射花样^[21], 如图 8所示。需要说明的是, 这里的计算假定了传播介质是全空间各向同性, 而震源区是各向异性介质。图 8中的断层面走向30°, 断层面倾角40°, 断层面滑动角为20°, 从上到下依次为各向同性介质、干燥HTI各向异性介质和饱和HTI各向异性介质, 从左到右张裂角分别为0, 10°, 50°, 90°。

4 结束语

本文分析了含裂隙干燥和饱和HTI介质的弹性参数以及各向异性HTI介质中的剪张源矩张量与微地震矩张量, 比较了震源区为各向同性介质、干燥和饱水HTI各向异性介质情形下的剪张源震源机制与地震矩张量沙滩球, 研究了震源区HTI各向异性对双力偶分量(DC)、补偿线性偶极子分量(CLVD)和各向同性分量(ISO)的影响以及各向异性HTI介质中剪张源的辐射花样。研究结果表明, 各向异性介质中的矩张量参数依赖于各向异性的类型和强度以及断层面的方向等参数。研究各向异性矩张量能够提供震源区的各向异性信息, 有利于了解地震如何破裂和扩展, 具有重要的价值和应用前景。