碳酸盐岩的矿物组成、结构以及成岩演化过程非常复杂^[1]。碳酸盐岩储集空间类型及其孔缝尺度的多样性导致该类储层具有强烈的非均质特征, 其不同类型岩相的波阻抗等地球物理参数呈现交叉叠置的现象。因此直接进行目标层系的地震波阻抗反演难以有效预测和描述储层特征^[2]。长期以来, 碳酸盐岩储层的地震精细预测一直困扰着地球物理工作者。

碳酸盐岩的发育及其储层特征总体上受沉积环境、成岩作用及其改造作用的影响, 尤其是不同沉积相带的储集性能和岩石物理特征均存在差异^[3-4], 因此, 刻画沉积微相及其展布对碳酸盐岩储层的预测具有重要的指示意义^[5-6]。为进一步提高复杂储层的地震预测精度, 充分考虑沉积相带的影响, 在地震相模式约束下将非线性随机模拟技术与地震反演相结合, 使地震反演储层预测结果更符合实际情况。HUANG等^[7]、KORNEEV等^[8]、ZHANG等^[9]和ZHOU等^[10]基于非线性理论及统计学方法, 探索并提出了高精度非线性混沌反演算法, 有效降低了地震反演的多解性, 从本质上提高了地震预测的效果。慎国强等^[11]通过叠前叠后的联合反演进一步提高了地震预测精度。黄捍东等^[12]、张志伟等^[13]、郭恺等^[14]、张显文等^[15]、杨子川等^[16]和李斌等^[17]针对塔河、阿姆河等不同地区的碳酸盐岩储层, 在地震地质特征分析的基础上, 建立岩相与地震反射特征的关系, 而后采用相带约束的方法在等时格架的基础上开展地震属性特征分析及相控反演研究, 改进了地震反演方法, 从而进一步提高了储层预测的精度。关达等^[18]针对大牛地气田D井区储层非均质性强的特点, 采用相控非线性地震反演技术, 一定程度上克服了运算结果的多解性问题, 该技术在华北油田复杂储层地震预测中发挥了重要作用, 取得了良好的应用效果。

本文研究区的盐下碳酸盐岩储层埋深大、非均质性强, 储层与非储层波阻抗值域范围严重重叠。为提高复杂碳酸盐岩储层的地震预测精度, 首先在岩石物理分析及地震相特征研究的基础上, 根据实际地震资料建立包含目标层段的宏观相控模型, 然后在测井资料的约束下采用非线性随机算法逐道、逐次地进行外推反演, 再进一步对反演数据体进行地震精细解释, 最终结合神经网络技术实现储层的定量预测。

1 方法原理

针对复杂的碳酸盐岩储集体, 采用地震与沉积相相结合的技术思路, 将地质研究成果与地震反演有机结合, 使地震反演结果最大限度地排除边界效应, 减少储层地震预测的多解性。地震相控逐道外推反演方法是从井旁道出发进行反演, 在精确的层序地层格架控制下实现反演, 并将每次得到的反演结果作为下一道反演的初始模型, 以此类推。

由于地震道与波阻抗的关系是非线性的, 所以基于地震道的逐道外推反演也可称为非线性随机反演。实际操作时可依据井的波阻抗和地震相划分结果建立地下波阻抗模型, 然后利用井点波阻抗和井旁地震道求出控制参数, 即可进行井约束、相控非线性随机反演。在实际资料反演过程中, 该反演方法可以有效地提高地震资料的垂向分辨率, 并充分考虑目标地层层序的地质特征, 使反演结果更符合实际情况。由此, 基于地震道非线性最优化反演思想, 假设岩石密度为常数, 将波阻抗反演转化为速度反演, 建立地震道与波阻抗关系的目标函数为:

$ f(V) = \sum\limits_{i = 0}^{n - 1} {{{(S_i^\Delta - {D_i})}^2}} \to \min $

(1)

式中:V为速度; S_i^Δ为模型响应, 即速度预测结果对应的合成地震记录; D_i为实际地震记录, i为地震记录的采样点序号, 范围为[0, n-1]。

利用Taylor公式将(S_i^Δ－D_i)在初始模型响应S_i处展开, 为方便求解将二次项以上的高次项省略, 可得:

$ \begin{array}{*{20}{c}} {S_i^\Delta - {D_i} = {S_i} - {D_i} + \sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{\partial {S_i}}}{{\partial {V_k}}} + }\\ {\frac{1}{2}\Delta V_k^2(\sum\limits_{k = 0}^{n - 1} {\frac{{{\partial ^2}{S_i}}}{{\partial {V_k}\partial {V_j}}}} ) + \cdots } \end{array} $

(2)

式中:S_i代表速度初始模型对应的合成地震记录; ΔV_k, ΔV_j代表不同模型参数的摄动量; k, j分别代表模型数和参变量。

对(2)式求导, 可得:

$ \begin{array}{*{20}{c}} {\frac{{\partial S_i^\Delta }}{{\partial \Delta {V_j}}} = \frac{{\partial {S_i}}}{{\partial {V_j}}} + \sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{{\partial ^2}{S_i}}}{{\partial {V_j}\partial {V_k}}}}\\ {(i = 0, 1, \cdots , n - 1\quad j = 0, 1, \cdots , n - 1)} \end{array} $

(3)

对(1)式的ΔV求导并展开, 可得:

$ \sum\limits_{i = 0}^{n - 1} {\frac{{\partial {{(S_i^\Delta - {D_i})}^2}}}{{\partial \Delta {V_j}}}} = 2\sum\limits_{i = 0}^{n - 1} {(S_i^\Delta - {D_i})} \frac{{\partial S_i^\Delta }}{{\partial \Delta {V_j}}} = 0 $

(4)

将(2)式和(3)式代入(4)式, 则有:

$ \begin{array}{*{20}{l}} {2\sum\limits_{i = 0}^{n - 1} \{ [{S_i} - {D_i} + \sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{\partial {S_i}}}{{\partial {V_k}}} + \frac{1}{2}\Delta V_k^2 \cdot }\\ {\quad (\sum\limits_{k = 0}^{n - 1} {\frac{{{\partial ^2}{S_i}}}{{\partial {V_k}\partial {V_j}}}} )](\frac{{\partial {S_i}}}{{\partial {V_j}}} + \sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{{\partial ^2}{S_i}}}{{\partial {V_j}\partial {V_k}}})\} = 0} \end{array} $

(5)

将(5)式两边除以2, 再展开, 可得:

$ \begin{array}{l} \mathop \sum \limits_{i = 0}^{n - 1} \left[ {({S_i} - {D_i})\frac{{\partial {S_i}}}{{\partial {V_j}}} + \frac{{\partial {S_i}}}{{\partial {V_j}}}\sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{\partial {S_i}}}{{\partial {V_k}}} + } \right.\\ \;\;\;\left. {({S_i} - {D_i})\sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{{\partial ^2}{S_i}}}{{\partial {V_j}\partial {V_k}}}} \right] + \mathop \sum \limits_{i = 0}^{n - 1} \left[ {\sum\limits_{k = 0}^{n - 1} \Delta {V_k} \cdot } \right.\\ \;\;\;\frac{{\partial {S_i}}}{{\partial {V_k}}} \cdot \sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{{\partial ^2}{S_i}}}{{\partial {V_j}\partial {V_k}}} + \frac{1}{2}\Delta V_k^2(\sum\limits_{k = 0}^{n - 1} {\frac{{{\partial ^2}{S_i}}}{{\partial {V_k}\partial {V_j}}}} ) \cdot \\ \;\;\;\left. {(\frac{{\partial {S_i}}}{{\partial {V_j}}} + \sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{{\partial ^2}{S_i}}}{{\partial {V_j}\partial {V_k}}})} \right] = 0 \end{array} $

(6)

在(6) 式中, 省略部分高阶极小量, 即省去左边第2项表达式, 得到(6)式的简化形式如下:

$ \begin{array}{l} \mathop \sum \limits_{i = 0}^{n - 1} \left[ {({S_i} - {D_i})\frac{{\partial {S_i}}}{{\partial {V_j}}} + \frac{{\partial {S_i}}}{{\partial {V_j}}}\sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{\partial {S_i}}}{{\partial {V_k}}} + } \right.\\ \;\;\;\;\;\;\left. {({S_i} - {D_i})\sum\limits_{k = 0}^{n - 1} \Delta {V_k}\frac{{{\partial ^2}{S_i}}}{{\partial {V_j}\partial {V_k}}}} \right] = 0 \end{array} $

(7)

为便于理解, 将(7)式用简单形式表示为:

$ A\Delta V + B\Delta V + C = 0 $

(8)

其中,

$ \left. \begin{array} { l } { A = \sum\limits_ { i = 0 } ^ { n - 1 } ( \frac { \partial S _ { i } } { \partial V _ { j } } \sum\limits_ { k = 0 } ^ { n - 1 } \frac { \partial S _ { i } } { \partial V _ { k } } ) } \\ { B = \sum\limits_ { i = 0 } ^ { n - 1 } [ ( S _ { i } - D _ { i } ) \sum\limits_ { k = 0 } ^ { n - 1 } \frac { \partial ^ { 2 } S _ { i } } { \partial V _ { j } \partial V _ { k } } ] } \\ { C = \sum\limits_ { i = 0 } ^ { n - 1 } [ ( S _ { i } - D _ { i } ) \frac { \partial S _ { i } } { \partial V _ { j } } ] } \end{array} \right. $

(9)

利用(8)式求取模型参数摄动量ΔV时, 常采用矩阵求逆的方法, 该方法容易因矩阵奇异而造成无解。为此, 可将矩阵求逆蜕变为一元一次方程求解以减少反演的多解性, 增强其稳定性。

将求得的ΔV代入(10)式, 多次迭代后可获得最终的反演速度V。

$ {V^{m + 1}} = {V^m} + \Delta {V^m} $

(10)

式中:m为迭代次数。

一般而言, 理想的地震反演结果需要地震相控模型的良好匹配和复杂运算的多次迭代。基于层序特征的地震相模型控制了地震反演的低频趋势, 为地震反演提供了必要条件; 非线性递推反演是实现地震反演的技术关键。

2 技术应用 2.1 研究区地质概况

南大西洋S盆地的L油田位于海域深水区, 含油气主要目标层为盐下白垩纪断凹期的BVE100、BVE200和BVE300层段。钻井揭示, 构成该套沉积地层的岩性主要包括微生物灰岩、白云岩、泥质灰岩、泥岩以及盐岩、膏岩等。根据岩心观测结果, 结合测井资料将储层分为3类:一类储层主要发育于礁滩相, 二类储层多分布于滩缘相, 三类储层基本为滩间洼地和半深湖-深湖相沉积岩。其中, 礁滩相主要包括叠层石、颗粒灰岩或球粒灰岩, 储层物性好, 孔隙度大于10%;滩缘相主要发育微晶灰岩和层纹岩等, 储层物性中等, 孔隙度为5%~10%;滩间洼地和半深湖-深湖相主要发育泥灰岩及颗粒灰岩, 储层物性相对较差, 孔隙度小于5%。因此, 有效预测此类盐下微生物灰岩优质储层的分布对深水区盐下油气勘探开发意义重大。

2.2 实际应用及效果分析 2.2.1 层序划分及地震相特征

L油田地震反射资料受上覆巨厚盐岩的影响, 目的层地震反射波视主频约为20 Hz。结合区域认识和钻井资料, 将该区盐下目标层段划分为4个三级层序:BVE组、ITP组、PIC组、CAM组, 并将重点目标BVE组进一步细化为3个四级层序:BVE100、BVE200、BVE300(图 1)。

分层系展开地震属性优化分析, 研究层系单元内部属性变化特征, 进一步综合构造、沉积以及古地貌等特征, 认为L油田研究区发育微生物礁、滩相复合体以及滩间洼地和半深湖相沉积。其中, 微生物礁的主体发育藻叠层石微生物灰岩, 局部与颗粒灰岩互层, 礁体周缘主要发育鲕粒灰岩和颗粒灰岩等, 成层性较好, 储层物性最好; 微生物滩相复合体因水动力变弱, 发育一些相对高速的致密灰岩和层状灰岩等, 储层物性中等; 滩间洼地和半深湖相沉积发育灰岩和泥质灰岩, 储层物性相对较差。

2.2.2 岩石物理特征

对研究区的4口井进行分层段岩石物理多参数交会分析, 主要将密度和纵波速度分别与纵横波速度比(v_P/v_S)、泊松比(μ)和拉梅系数乘以密度(ρ)交会, 分析岩性与弹性参数关系, 统计岩石物理特征。如图 2所示, 以BVE100层段为例, 其中一类和二类储层的纵横波速度比、泊松比和拉梅系数与密度的乘积均较低, 一类和二类储层的纵横波速度比、泊松比和拉梅系数与密度的乘积存在一定范围的重叠; 三类储层的纵横波速度比、泊松比和拉梅系数与密度的乘积均较高。BVE100层段不同类型储层的弹性参数统计结果见表 1。

2.2.3 相控地震反演

研究区三维地震资料面积240 km², 选定其中的7口井参与分析和运算, 确定地震反演的3个目标层段, 即4个约束层段界面(BVE100、BVE200、BVE300及ITP层段顶界面), 以地层层序界面和前期获得的地震相模型为约束展开地震相控非线性反演。钻井揭示, BVE100层段发育叠层石灰岩及颗粒灰岩组合, 与上覆厚层盐岩相比, 具有明显的低波阻抗特征; BVE200层段是大套的泥灰岩夹薄层灰岩, 与上、下围岩的地震反射界面清晰; BVE300层段与下伏地层存在明显的角度不整合特征, 主要发育大套泥灰岩组合, 其波阻抗值通常较大。相控地震反演过程如下:首先通过井震综合标定, 建立基于层序特征的地震相控低频模型, 并将其作为反演的基础条件; 然后以模型参数作为区域化变量, 采用上述非线性逐道外推方法展开随机模拟; 再调整变差函数以确立模型参数点之间的统计相关关系, 基于目标函数最小化原则将模型参数变量作为符合高斯分布的随机变量; 最后通过高精度随机模拟实现约束反演结果与地震道的最佳匹配, 进而获得研究区的地震反演速度、密度等参数。图 3为地震相控非线性反演的波阻抗与反演得到的地震剖面的叠合显示(红线分别为4个约束层段顶界面)。反演得到的波阻抗不仅与地震波组的振幅、频率、相位的细节变化高度一致, 而且反演得到的地震剖面中地震反射的中、弱信号更加清晰, 目标层内幕细节更加突出和丰富。

2.2.4 应用效果 2.2.4.1 岩性岩相和微生物礁特征

图 4为过8-LL-9井与过1-R-682A井的地震反演波阻抗剖面与井点处目标层段岩性岩相特征, 井震特征具有良好的一致性。图 4a中的灰岩储层和致密的泥灰岩、图 4b顶部的高速石膏和下部较高速度的白云岩波阻抗均存在明显差异, 波阻抗反演结果进一步提高了目标层段的纵、横向分辨率。由图 5a可知, 经地震相控非线性反演后, 微生物礁体显示得更清晰, 但受地震分辨率的限制, 从图 5b中无法辨识微生物礁体的细节, 只能看到一些较弱的内部反射特征。图 5a不仅更清晰地反映了微生物礁体的叠瓦状轮廓, 也进一步揭示了其内幕结构特征。

2.2.4.2 储层厚度分布预测

基于地震相控非线性反演结果, 利用纵波速度、波阻抗等参数与井点储层厚度进行多参数拟合分析, 最终预测出研究区储层厚度的分布情况, 进而编制出目标层段一类和二类储层以及三类储层的预测厚度平面展布(图 6)。从图 6a可以看出, BVE100层段一类和二类优质储层集中发育在两个构造隆起区, 最厚的储层发育在南部, 厚度达70 m, 其次是北部的东西相邻井区, 储层最大厚度分别为50 m和25 m, 平均厚度约为20 m。从图 6b可以看出, BVE100层段局部发育三类储层, 厚度较小, 主要分布于研究区的中部, 呈条带状南北向局限分布, 储层平均厚度约10 m。

对研究区14口井依次统计3个层段的地震反演预测厚度与测井解释厚度吻合情况, 将二者相对误差大于20%视为不吻合, 统计结果如表 2所示。储层厚度预测的总吻合率为90.9%, 表明地震相控非线性反演技术能够有效实现复杂碳酸盐岩储层的精细预测。

3 结论认识

1) 相控地震反演技术的核心思想是将地震相时空分布模型与地震反演运算有机结合, 使地震反演最大限度地排除边界效应。理想的地震反演需要相控模型的良好匹配和复杂运算的多次迭代。其中的基于层序特征的地震相模型控制了反演的低频趋势, 为地震反演提供了必要条件。

2) 相控低频模型是实现复杂碳酸盐岩储层精细预测的重要基础, 基于相界面和测井资料建立的低频速度模型理想地控制了岩性岩相的宏观变化, 从根本上减少了地震反演的多解性。

3) 实际资料应用统计结果表明, 后验井符合率大于90%, 储层厚度与测井资料的相对误差小于20%。该技术具有运算效率高、反演结果随机性小的特点, 可以更好地体现相控的思想, 并有效实现复杂储层的定量预测。