摘要 针对目前用于建立二维速度结构模型的 RAYINVR 软件对四分量海底地震仪(OBS)记录的转换横波无法自动反演, 建模效率低的问题, 采用 MATLAB 遗传反演算法, 对 RAYINVR 软件进行改进, 实现对横波速度结构模型各层和层中区块泊松比的自动同步反演, 可以为天然气水合物勘探调查提供泊松比和杨氏模量等重要弹性参数。利用采自 Slipstream 海底滑坡的 OBS 数据, 通过同步反演方法, 获得该地区较为精细的纵波和横波速度结构模型, 并可与附近 U1326 钻孔的测井数据进行对比, 验证了横波速度结构同步反演建模方法的有效性。模型揭示了两个泊松比差异大的结构面: 代表水合物稳定底界的似海底反射界面(BSR)(海底之下230±10m)以及浅部异常高速体(有可能是高饱和度水合物富集的砂体)的底界(海底之下 75～100m), 后者与滑脱面大致重合, 指示水合物与Slipstream海底滑坡的形成有关。
关键词 同步反演; 横波速度结构; Slipstream海底滑坡; 天然气水合物; 射线追踪
地震勘探是海域天然气水合物调查的重要途径[1‒2], 不但可以对水合物储层进行成像刻画, 还可以建立速度结构模型, 从而有效地评估水合物资源量, 对后续开采意义重大。常规的海面拖缆地震勘探仅能产生和接收单一的纵波信号, 对沉积物孔隙流体中的气体较为敏感, 少量的游离气就可能造成大幅度的纵波能量衰减和速度降低[3], 影响水合物地层的成像效果和饱和度估算准确度。这种现象在作为水合物稳定区底界标志的地震似海底反射层(bottom simulating reflector, BSR)上下(表现为杂乱强反射)或携带大量游离气向上运输的冷泉通道(气烟囱)内部(表现为条带状空白)尤为严重[4]。
横波主要通过沉积物骨架传播, 几乎不受孔隙流体的影响[5]。因此, 在海域水合物重点勘探区, 通常使用四分量海底地震仪(Ocean Bottom Seismo-meter, OBS)来记录由纵波在地层界面上转换而来的横波信号[6‒7], 从而建立沉积地层的纵横波速度结构模型, 以便获取更加丰富、更加准确的水合物储层信息。
加拿大西海岸温哥华岛外的北卡斯卡底(Casca-dia)陆缘增生楔主要由板块俯冲过程中刮削下来的洋盆沉积物堆积而成[8‒9], 蕴藏着丰富的天然气水合物资源[10‒11], 20 世纪 80 年代以来一直是水合物研究的重点区域, 不但具有丰富的单道地震(single channel seismic, SCS)、多道地震和 OBS 等地球物理勘探数据[12‒13], 而且有专门针对水合物的大洋钻探计划(Ocean Drilling Program, ODP) 146 航次[14]和综合大洋钻探计划(Integrated Ocean Drilling Program, IODP)X311 航次[15‒17]的钻井数据。特别之处是, 在该增生楔前缘的条状变形脊上发育一系列海底滑塌构造。Lopez 等[18]和 Yelisetti 等[19]分别运用 OBS 数据和 RAYINVR 软件[20], 研究其中规模较大的 Orca和 Slipstream 海底滑塌的地震速度结构, 认为其形成与水合物富集层的位置有关。然而, 受限于转换横波处理的复杂性和 RAYINVR 软件的缺陷, 他们只建立滑塌构造的纵波速度模型, 缺乏横波速度信息, 无法提供滑塌构造力学分析所需的杨氏模量和泊松比等弹性参数。
广泛使用的 RAYINVR 软件是用射线追踪方法对二维速度模型进行走时拟合, 其内置的阻尼最小二乘(damped least squares)反演算法能够方便地建立纯纵(横)波速度结构模型。但是, 对于 OBS 记录的转换横波反射事件, 该软件只能通过人工改变模型的泊松比参数, 将纵波速度转换为相应的横波速度, 并且, 一旦存在横向速度结构不均的现象, 还需要将该结构层进一步划分区块, 模型会变得更加复杂。因此, 横波模型需要花费大量的时间进行反复的手动测试, 效率过低, 严重妨碍 OBS 横波信号的处理。
针对上述问题, 本研究在 RAYINVR 软件的基础上, 依托 MATLAB 语言建立 OBS 横波模型的遗传反演方法, 实现各结构层和层中区块参数(泊松比)的同步反演。利用该方法建立的较精细的 Slips-tream 海底滑坡纵、横波速度结构模型, 指示浅部高速异常体(海底之下 75～100m)和深部 BSR (海底之下 230±10m)两个软弱面为潜在的滑脱面位置。
对 OBS 等广角地震数据, 通常采用射线路径追踪的算法进行速度结构层析反演, 采用正演试错拟合的求解方案对模型参数进行迭代[21‒22], 直至模型计算的走时与真实记录的地震信号走时之间误差足够小。对于二维速度结构模型, 目前普遍使用 Zelt等[20,23]编写的 RAYINVR 软件, 其模型构建采用不规则四边形网格, 通过灵活地设置速度以及界面节点的个数和位置, 使得模型层界面与地震事件识别情况一致, 能够在很大程度上避免过参数化引起的反演假象以及欠参数化引起的细结构丢失问题。模型可以分为若干个不规则四边形, 预先定义每个四边形的顶点速度, 然后采用双线性内插法, 由顶点速度得到四边形内部的速度, 因此在四边形内部的射线追踪存在解析解。该软件可以对折射波、反射波、界面滑行波、转换波和层间反射波等进行射线追踪计算, 具有运算速度快和算法稳定等优点。
当没有转换波时, RAYINVR 软件可以使用自带的阻尼最小二乘反演工具(DMPLSTSQR 程序), 从上往下逐层对模型节点速度和(或)层面深度参数进行自动的快速的迭代反演, 每一层的射线走时拟合流程见图 1。以建立纵波速度结构模型为例, 首先用 RAYINVR 软件正演计算某层的射线路径和走时, 需要输入 v.in, tx.in 和 r.in 这 3 个文件。其中, 文件 v.in 包含模型每一层的深度和速度信息, 文件tx.in 是从地震成像数据拾取的真实地震事件走时, 文件r.in则是关于射线追踪的算法控制以及结果画图的参数设置。RAYINVR 正演会输出射线追踪走时与实际地震事件走时残差的均方根(root mean squares, RMS)和标准卡方值(normalized χ2), 后者的计算公式为
其中, N 为实际拾取的地震走时数量, Ti 为第 i 个实际拾取的地震走时, ti 为相应的射线追踪得到的地震走时, σi 为第 i 个地震走时的拾取误差。在得到标准卡方值后, 通过人机交互判断 χ2 是否满足要求: 当 χ2≤1 时, 说明走时残差达到地震事件拾取误差水平, 可以结束反演; 如果 χ2 显著大于 1, 则将该值输入 DMPLSTSQR 程序, 同时输入反演参数控制文件d.in (在该文件中选择需要反演的目标界面), 通过偏微分方法求出模型参数调整向量, 据此更新模型文件 v.in 中的速度和(或)深度参数。迭代上述正演和反演过程, 直到获得模型最优解(保存在最后一步的文件 v.in 中)。
可见, DMPLSTSQR 程序是通过改变文件 v.in 中的速度和地层界面深度参数来实现迭代反演的。然而, 当对转换横波进行射线追踪时, 纵波转换为横波的反射界面参数(如转换发生的层号和泊松比等)均存储在文件 r.in 中, 无法由 DMPLSTSQR 程序自动更新, 只能手动修改。因此, 横波速度模型的迭代反演过程非常耗时, 并且, 当研究区的地质构造比较复杂时, 需要对目标层进一步细化分块, 以至很难手动将模型参数调节到令人满意的程度。
图1 RAYINVR软件迭代反演速度模型流程 (修改自文献[20])
Fig. 1 Flow chart of iterative velocity model inversion by RAYINVR software (modified from Ref. [20])
根据震动形式, 可将地震波分为体波和面波, 其中体波又分为纵波和横波。当纵波到达一个反射界面时, 其反射波中不仅包含纵波, 还包含转换而来的横波(图 2)。由于海面震源激发的地震波在海水中只能传播纵波, 因此放置在海底的 OBS 记录的主要横波事件(能量最强, 可以识别)是由纵波在其反射界面上转换产生。因此, RAYINVR 对 OBS 横波的射线追踪, 不是直接确定横波速度, 而是在给定的纵波速度结构模型(地层界面深度和纵波速度梯度保持固定)基础上, 寻找最优的地层泊松比。在通过 RAYINVR 软件得到模型的泊松比 v 之后, 横波速度VS则由对纵波速度VP的变换得到:
。如果研究区域的地质构造较为简单, 可以通过手动调节的方法, 给模型各层设置泊松比参数。但是, 当地质构造比较复杂时, 地层的泊松比不仅在纵向上可能存在差异, 在横向上也可能存在较大的差异, 需要在层内进一步划分区块。此时, 用人工调节参数的方法就会事倍功半, 因此需要研发一套能够让程序实现自动同步反演而得到最优解的建模方法。
θ1为P波在反射界面的入射角, θ2为P波的折射角, ϕ1为转换S波的反射角, ϕ2为转换S波的折射角
图2 纵波在弹性界面发生反射和折射示意图
Fig. 2 Reflection and refraction diagram of P wave incident on an elastic interface
遗传算法通过模拟自然进化过程, 对值域空间进行搜索, 得到全局最优解, 是一种可用于复杂系统优化的算法。本文将 RAYINVR 软件的射线追踪法与 MATLAB 软件的遗传算法相结合, 用 MATLAB语言编写一套能够实现 OBS 横波速度结构模型在各层以及层间不同区块同步反演的方法, 从而节省人工修改参数的时间成本, 实现复杂地质构造环境下反演横波速度结构的自动化。该同步反演方法的基本流程见图 3。
1)给定模型各层和层中区块的初始泊松比组合(保存在文件 r.in 中)。
2)设定遗传算法的最大代数 N 和每代的种群数量 M, 启动自动反演。
3)RAYINVR 软件读入射线追踪控制文件 r.in、拾取的地震走时文件 tx.in 以及速度和层界面深度控制文件 v.in, 通过正演射线路径拟合得到 χ2。
图3 OBS横波速度模型同步反演流程
Fig. 3 Flow chart of S-wave velocity modeling for OBS records by synchronous inversion method
4)遗传算法会判断此时是否达到预设条件: 如果达到预设的迭代次数 N, 或者 χ2≤1 (意味着此时的泊松比为全局最优解), 反演终止; 如果没有达到预设的迭代次数 N, 并且 χ2>1, 那么程序会生成一组新的泊松比组合, 并自动传给 RAYINVR 软件, 进行新一轮正演。
5)重复上述迭代过程, 直到满足预设条件, 得到的最新泊松比组合将保存在文件 r.in 中。
在卡斯卡底活动陆缘, 胡安德富卡(Juan de Fuca)板块以约 46mm/a 的速率俯冲到北美板块之下, 来自洋盆的大量沉积物在变形前缘位置被刮削下来, 堆积形成增生楔[8]。这些增生体沿着一些巨大的深达洋壳的逆冲断层发生挤压和缩短, 并相互重叠, 从而在陆坡底部(变形前缘)形成一系列与板块边界近平行的长条状背斜山脊[19], 本文的研究对象 Slip-stream 山脊就是其中之一(图 4)。卡斯卡底活动陆缘一直是天然气水合物研究的聚焦区域, 因为大量地震调查都观测到从变形前缘到陆坡中部的深水沉积物中广泛存在作为水合物稳定区底界地震标志的BSR, 海底之下的深度约为 200～300m, 显示大量天然气水合物的富集[5,24‒26,28], 也得到 ODP146 和 IODP X311 大洋钻探结果的证实。
近年来, 在高分辨率多波束海底地貌的调查中发现, 变形前缘的条状背斜山脊上经常发育滑塌构造。例如, IODP X311 U1326 站位所在的山脊发育一个规模较大的 Orca 海底滑塌(图 4), 其滑脱面深度恰好抵达 BSR, 因此推测该滑塌构造的形成与水合物层有关[18]。作为对比, 其相邻的 Slipstream 滑塌(U1326 站位东南方向约 15km)规模更大, 推测其形成也与水合物有关, 但滑脱面的深度却只有 100m左右[19]。
根据美国华盛顿大学 TN175 航次采集的 EM300多波束测深数据[19], Slipstream 海底滑坡体的长度约为 3.3km, 宽度约为 2.5km, 滑坡角度较大(约为25°), 出露于海底的陡壁高度可达近 100m。该滑坡体上发育一系列正断层, 且部分断层面出露于海底, 因此滑坡体在地貌上类似多个排列不整齐的梯形箱状高地。较大的滑坡块体主要分布在距离山脚 1～2km 的海底平原沉积层之上, 一些粒度较小的滑坡体沉积物则外流到更远的地区(图 4)。本研究采用的 SCS 和 OBS 数据来自 SeaJade (Seafloor Earthquake Array – Japan Canada Cascadia Experiment)项目。
修改自文献[27]。红色虚线框为 Slipstream 海底滑坡区域; 红色三角形为 OBS 测点 30, 31, 33 和 34 的位置; 4 号地震测线平行地穿过Slipstream 山脊, 其中绿线为 OBS 地震航次, 黄线为 SCS 航次; 红色五角星为 IODP X311 航次的U1326 站位; 多波束海底地貌图上叠加的等深线间隔为 25 m
图4 Slipstream海底滑坡区构造背景和多波束海底地貌
Fig. 4 Tectonic background and multibeam bathymetry of Slipstream submarine slide
由于 OBS 记录的横波是由纵波在地层界面上转换而来, 因此需要建立 Slipstream 海底滑坡区域的纵波速度结构模型。
传统的海面拖缆 SCS 数据能够很好地针对海底沉积层的结构特点进行成像, 但受限于炮点与接收器之间的偏移距(source-receiver offset)在深海中太小, 无法得到准确的速度和深度信息。布设于海底的 OBS 采集的是广角(大偏移距)数据, 能够很好地对海底沉积层的速度信息进行分析, 但是对其结构信息的表征则有所欠缺。因此, 在构建海底沉积层纵波速度结构模型时, 我们同时采用 OBS 广角数据和 SCS 数据来实现对地层速度和深度的同步约束。需要注意的是, RAYINVR 软件最初是用于大陆勘探的, 设计的是单震源对应多个接收器, 而海洋的地震勘探则是单个 OBS 记录多个炮点数据, 因此在建立模型时, 需要将 OBS 作为射线发射点, 将炮点作为接收位置。
前人在卡斯卡底增生楔区域创建的纵波模型一般为 5～7 层[5,19,25], 用于滑塌区力学机制分析时精细度不够。我们通过仔细地辨识 OBS 和 SCS 的反射和界面滑行波事件(图 5), 在 SCS 数据中拾取 10 个纵波反射事件(line4A_r1～line4A_r10), 在 OBS 数据中拾取 7 个纵波反射事件(RL0～RL6)和 5 个纵波临界折射(界面滑行波)事件(RF1～RF5)。
(a)在 4 号地震测线的 SCS 数据中拾取的共 10 个纵波反射事件(line4A_r1～line4A_r10), 其中 r1 对应海底, r2 和 r3 分别对应浅部高速层的顶、底界面, r5 对应BSR; (b)在 OBS 测点 33 的数据中拾取的广角地震反射(RL)和临界折射事件(RF), 其中 RL0, RL1 和RL2 分别对应海底、高速层顶部和 BSR, RF1 对应高速层顶界面, RF2～RF5 为 BSR 之下的深部临界折射事件, 下同。折合走时=走时−偏移距/折合速度, 本文折合速度为 4000 m/s
图5 Slipstream滑坡区拾取的纵波反射和临界折射地震事件
Fig. 5 Reflection and critical refraction (head wave) events of P wave for Slipstream slide
综合 SCS 和 OBS 数据的拾取事件, 最终将 Slip-stream 海底滑坡沿 4 号测线方向的山脊剖面分为 14层(不包含海水)。反演得到的最终射线路径追踪结果见图 6(a), 拟合得到的走时与实际拾取的走时相符(图 6(b)), 模型整体的 RMS 走时残差为 15ms, 小于拾取误差(20ms), χ2=0.901, 达到预期要求。
精细化的纵波速度结构模型实现对最小厚度为25m 左右薄层的识别(图 7(a)), 其中 L1 与海水直接接触, 平均速度约为 1.51km/s, 代表非常松散的海底沉积物; L2 的顶部界面在海底之下约 75m, 厚度约为 25m, 纵波速度突然增加至约 2.0km/s, 推测可能是有高饱和度天然气水合物富集的砂体; L3 ～L5, 纵波速度从 1.7km/s 缓慢增加至 2.1km/s, 代表无(或含少量)水合物沉积物的压实趋势; L5 的底部界面对应 BSR, 深度约在海底之下 230m; BSR 之下的L6 纵波速度略有削减, 代表游离气体的存在; 之后, 纵波速度的增加再次向压实趋势线回归。
(a)纵波速度结构模型拟合的射线路径, 4个圆圈为OBS测点位置及编号(下同); (b)纵波速度结构模型的走时拟合结果
图6 OBS事件在Slipstream滑坡区的纵波射线追踪结果
Fig. 6 Ray tracing and travel-time fit results of OBS P-wave events for Slipstream slide
当沉积物孔隙中含游离气体时, 不但会影响SCS 反射数据和 OBS 垂直分量的成像数据, 衰减严重时还会出现地震反射空白, 大幅度地降低纵波速度, 影响水合物饱和度的计算。然而, 横波传播不受孔隙流体影响, 可以反映含气区域(如 BSR 之下的游离气层)的内部结构。并且, 横波速度在 BSR之上的大幅度增加可能显示固体天然气水合物的存在。此外, 由于 OBS 记录的横波是由纵波转换而来, 其频率与纵波相同, 但速度更慢(横波在纯天然气水合物中的速度约为 1890m/s, 在不含水合物的浅层海底沉积物中的速度为 100～600m/s)[29], 因此波长更短, 从而可以解析比纵波更精细的沉积结构信息。
在上述纵波速度结构模型的基础上, 保持各层顶、底界面的深度节点和纵波速度节点不变, 通过MATLAB 的遗传算法对横波模型的泊松比进行同步反演。模型初始的参考泊松比数值(图 8(c))通过15 km 之外的 IODP U1326D 钻孔的纵、横波测井结果计算得到:
(a)纵波速度结构模型, L2为高速层(纵波平均速度约为2.0 km/s), L5 底部界面为 BSR; (b)横波速度结构模型
图7 Slipstream滑坡区纵横波速度结构模型
Fig. 7 P- and S-wave velocity models of Slipstream slide
需要指出的是, 通过非线性的遗传算法得到的满足条件的泊松比组合是全局最优解。作为对比, RAYINVR 的 DMPLSTSQR 程序是线性优化, 在初始参数与模型最终值偏差较大或地质构造复杂时, 容易掉进局部解的陷阱, 造成反演结果不佳, 需要重新手动调节初始参数再尝试。图 9 展示转换横波的射线追踪路径和走时拟合结果, 得到的 RMS 走时残差为 27ms, 小于横波事件的拾取误差(30ms), χ2=0.903, 达到预期要求。
将横波速度结构模型(图 7(b))与纵波速度结构模型(图 7(a))进行对比, 可见浅部的高速异常体 L2的横波速度增加幅度(14.3%)明显小于纵波(24.5%), 推测该层的高饱和度天然气水合物很可能以固体填隙物的形式富集, 对沉积物骨架的力学性能增加不明显(泊松比与其上接触海水的 L1 相近)。BSR 之上的 L5 横波速度相对于其下游离气层的增加幅度(14.0%)明显大于纵波(5.3%), 说明纵波速度受到游离气的干扰, 没有显示出含水合物层的高速特征。图 8(a)和(b)以 OBS 测点 34 所在位置的一维深度‒速度剖面形式对比上述纵、横波速度的变化情况, 图8(c)的泊松比则显示浅部高速层(L1)与 BSR 之下游离气层相比, 泊松比明显偏高, 指示这两个部位为力学软弱层。L2 的泊松比大于游离气层, 可能是Slipstream滑塌体沿着相对更软的高速层底界发育的原因。作为对比, Orca 滑塌体旁边的 U1326D 测井数据显示浅部的高速异常体出现在海底之下 60～ 90m 处, 而 Slipstream 海底滑坡区的高速层是在海底之下 75～100m, 这可能是 U1326D 与 Slipstream滑坡区之间约 15km 的空间距离引起的差异。此外, U1326D 的浅部高速异常体泊松比为低值, 指示力学强度相对较高, 并不是一个软弱面, 因此其滑脱面出现在泊松比值对比更强烈的 BSR 位置(应力容易在力学性质相差较大的分界部位集中), 显示与Slipstream 滑塌体不同的水合物控制作用。
图8 Slipstream海底滑坡区纵横波速度结构模型与U1326D测井数据[30]对比
Fig. 8 Comparison of inverted P- and S-wave velocity models of Slipstream slide with IODP X311 U1326D well loggings[30]
1)将 RAYINVR 软件射线追踪算法与 MATLAB遗传反演相结合, 实现 OBS 横波速度结构模型各层和层中区块泊松比的同步反演, 极大地提高了地质构造复杂区域的横波速度结构建模效率。
2)建立 Slipstream 海底滑坡体的精细纵、横波速度结构模型, 显示作为水合物稳定区底界地震标志的似海底反射层(BSR)位于海底之下 230～250m, 而浅部约 75～100m 的高速异常体很可能是高饱和度天然气水合物富集的砂体。模型结果可与 IODP X311 航次 U1326D 测井数据进行对比, 验证了同步反演速度结构建模的有效性。
3)浅部高速异常体的深度与滑坡面大致重合, 显示水合物造成的力学软弱面对滑坡构造的控制作用, 为后续建立三维力学模型提供了关键证据。
(a) OBS 横波速度结构模型拟合的射线路径, S1～S8代表转换横波反射事件; (b) OBS 横波速度结构模型的走时拟合结果
图9 Slipstream滑坡区转换横波的射线追踪结果
Fig. 9 Ray tracing and travel-time fit for S-wave velocity model of Slipstream slide
致谢 衷心感谢加拿大维多利亚大学 George Spence 教授提供 OBS, SCS 和多波束等数据, 特别感谢 IHS 全球公司提供地球科学综合解释软件包“The Kingdom Suite”。
参考文献
[1] Makogon Y F. Natural gas hydrates — a promising source of energy. Journal of Natural Gas Science and Engineering, 2010, 2(1): 49‒59
[2] Helgerud M B, Dvorkin J, Nur A, et al. Elastic-wave velocity in marine sediments with gas hydrates: effective medium modeling. Geophysical Research Letters, 1999, 26(13): 2021‒2024
[3] Yang J, Davies R J. Gravity-driven faults: migration pathways for recycling gas after the dissociation of marine gas hydrates. Marine Geology, 2013, 336: 1‒9
[4] Davies R J, Clarke A L. Methane recycling between hydrate and critically pressured stratigraphic traps, offshore Mauritania. Geology, 2010, 38(11): 963‒966
[5] Dash R, Spence G. P-wave and S-wave velocity struc-ture of northern Cascadia margin gas hydrates. Geo-physical Journal International, 2011, 187(3): 1363‒ 1377
[6] 张光学, 徐华宁, 刘学伟, 等. 三维地震与OBS联合勘探揭示的神狐海域含水合物地层声波速度特征. 地球物理学报, 2014, 57(4): 1169‒1176
[7] 刘伊克, 朱伟林, 米立军, 等. 南海深水多次波成像. 中国科学: 地球科学, 2015, 45(2): 152‒160
[8] Davis E E, Hyndman R D. Accretion and recent defor-mation of sediments along the northern Cascadia sub-duction zone. Geological Society of America Bulle-tin, 1989, 101(11): 1465‒1480
[9] Hyndman R D. The Lithoprobe corridor across the Vancouver Island continental margin: the structural and tectonic consequences of subduction. Canadian Journal of Earth Sciences, 1995, 32(10): 1777‒1802
[10] Mackay M E, Jarrard R D, Westbrook G K, et al. Origin of bottom-simulating reflectors: geophysical evidence from the Cascadia accretionary prism. Geo-logy, 1994, 22(5): 459
[11] Pohlman J W, Kaneko M, Heuer V B, et al. Methane sources and production in the northern Cascadia mar-gin gas hydrate system. Earth and Planetary Science Letters, 2009, 287(3/4): 504‒512
[12] Riedel M, Tréhu A, Spence G. Characterizing the thermal regime of cold vents at the northern Cascadia margin from bottom-simulating reflector distributions, heat-probe measurements and borehole temperature data. Marine Geophysical Researches, 2010, 31(1/2): 1‒16
[13] Hyndman R D, Riedel M, Spence G D. Hydrate studies of northern Cascadia margin off Vancouver Island: a reference source // 6th International Con-ference on Gas Hydrates (ICGH 2008). Vancouver, 2008, doi: 10.14288/1.0040999
[14] Westbrook G K, Carson B, Musgrave R J. Proceeding of Ocean Drilling Program, volumn 146 part 1: initial reports [R]. Texas: Ocean Drilling Program, Texas A&M University, 1994
[15] Torres M E, Tréhu A M, Cespedes N, et al. Methane hydrate formation in turbidite sediments of northern Cascadia, IODP Expedition 311. Earth and Planetary Science Letters, 2008, 271: 170‒180
[16] Riedel M, Collett T S, Malone M J. Gas hydrate drilling transect across northern Cascadia margin — IODP Expedition 311. Geological Society London Special Publications, 2009, 319(1): 11‒19
[17] Expedition 311 Scientists. Cascadia margin gas hyd-rates [R]. Washington, DC: Integrated Ocean Drilling Program Management International, 2005
[18] Lopez C, Spence G, Hyndman R, et al. Frontal ridge slope failure at the northern Cascadia margin: margin-normal fault and gas hydrate control. Geology, 2010, 38(11): 967‒970
[19] Yelisetti S, Spence G D, Riedel M. Role of gas hydrates in slope failure on frontal ridge of northern Cascadia margin. Geophysical Journal International, 2014, 199(1): 441‒458
[20] Zelt C A, Smith R B. Seismic traveltime inversion for 2-D crustal velocity structure. Geophysical Journal International, 1992, 108(1): 16‒34
[21] Catchings R D, Mooney W D. Crustal structure of the Columbia Plateau: evidence for continental rifting. Journal of Geophysical Research: Solid Earth, 1988, 93(B1): 459‒474
[22] Boland A V, Ellis R M. Velocity structure of the Kapuskasing Uplift, northern Ontario, from seismic refraction studies. Journal of Geophysical Research: Solid Earth, 1989, 94(B6): 7189‒7204
[23] Zelt C A, Ellis R M. Practical and efficient ray tracing in two-dimensional media for rapid traveltime and amplitude forward modeling. Canadian Journal of Exploration Geophysics, 1988, 24(1): 16‒31
[24] Hyndman R D, Spence G D. A seismic study of methane hydrate marine bottom simulating reflectors. Journal of Geophysical Research: Solid Earth, 1992, 97(B5): 6683‒6698
[25] Chen M. Northern Cascadia marine gas hydrate: con-straints from resistivity, velocity, and AVO [D]. Vic-toria: University of Victoria, 2006
[26] Lopez C. Seismic velocity structure associated with gas hydrate at the frontal ridge of Northern Cascadia Margin [D]. Victoria: University of Victoria, 2007
[27] 隆松伯, 何涛, 梁前勇, 等. 基于ANSYS Work-bench的高精度自动化三维地质建模方法——以天然气水合物相关的Slipstream海底滑坡为例. 北京大学学报(自然科学版), 2018, 54(5): 994‒1002
[28] Hyndman R D, Davis E E. A mechanism for the formation of methane hydrate and seafloor bottom-simulating reflectors by vertical fluid expulsion. Journal of Geophysical Research: Solid Earth, 1992, 97(B5): 7025‒7041
[29] Waite W F, Helgerud M B, Nur A, et al. Labora- tory measurements of compressional and shear wave speeds through methane hydrate. Annals of the New York Academy of Sciences, 1999, 912(1): 1003‒1010
[30] Riedel M, Collett T S, Malone M J. Expedition 311 summary [R]. Washington, DC: Integrated Ocean Drilling Program Management International, 2006
Model Inversion of Velocity Structure for Slipstream Submarine Slide Related to Gas Hydrate
Abstract The RAYINVR software, which is widely used in academia for 2-D velocity structure model from four-component OBS (ocean bottom seismometer) data, is unable to invert for the converted S-wave automatically, resulting in the low-efficiency of modeling process. Using MATLAB’s genetic algorithm, the RAYINVR software is improved and able to automatically and synchronously invert for Poisson’s ratios of each layer with all sub-blocks for the S-wave velocity structure model, and thus can provide Young’s modulus, Poisson’s ratio and other important mechanical information for gas hydrate survey. This method is applied to process the OBS data collected at the Slipstream submarine slide, and a fine P- and S-wave velocity structure model is obtained, which is comparable to the logging data of nearby borehole U1326. Therefore, the validity of the auto-synchronous inversion method is verified for the S-wave velocity structure modeling. The optimal velocity model reveals two structural interfaces with large Poisson’s ratio contrast. One is BSR (bottom simulating reflector) at 230±10 mbsf (meter beneath sea floor), which represents the bottom boundary of the gas hydrate stability zone, and the other is the basal boundary of a shallow abnormal high-speed body (possibly a sand body enriched with high saturation gas hydrate) at 75‒100 mbsf. The latter agrees roughly with the glide plane of Slipstream submarine slide, indicating that the hydrate is related to the formation of submarine landslide.
Key words synchronous inversion; S-wave velocity structure; Slipstream submarine slide; gas hydrate; ray tracing
doi: 10.13209/j.0479-8023.2021.005
国家自然科学基金(41676032)和中国地质调查局国家天然气水合物专项基金(DD20190234, DD20190218, DD20189320, DD20160217, HD-JJHT-20)资助
收稿日期: 2020–03–10;
修回日期: 2020–05–15