正演模拟技术共55页
利用正演模型研究VSP处理
利用正演模型研究VSP处理涂齐催1姜雨1杜本强2(1.中国海洋石油有限公司上海分公司研究院，上海200030)(2.冀东油田勘探开发研究院，唐山063004)摘要：VSP是一种井中地震观测技术，其独特的观测方式，决定了VSP处理也有别于常规地震，最明显的就是VSP资料必须经过波场分离，才能利用下行波或上行波信息对资料进行解释和反演。
除波场分离外，VSP还有许多其它独特的处理环节。
基于PC-VSP系统，进行了VSP模拟，在正演数据的基础上，研究了VSP处理，取得了较好的效果。
关键词：VSP处理；常规地震；波场分离；PC-VSP系统；VSP模拟Using forward model to study VSP processingTu Qicui1Jiang Yu1Du Benqiang2(1. CNOOC Shanghai Ltd., Shanghai 200030, China)(2. Exploration and Development l Research Institute, Jidong Oilfield, TangShan 063004, China)Abstract：VSP is a survey technology of borehole seismic. Because of its special survey, VSP processing is also different from surface seismic and the most different is that VSP data need to be processed with wavefield separation, after which, downgoing wave or upgoing wave can be used to do interpretation and further inversion. Besides, VSP has many other special processing steps. With PC-VSP system, VSP modeling is done and further research on VSP processing is done, based on forward data, and the final result of processing is very good.Key words：VSP processing；surface seismic；wavefield separation；PC-VSP system；VSP modeling──────────────── 1作者简介：涂齐催，1979年生，开发地震工程师，硕士研究生；现主要从事地震岩石VSP是一种井中地震观测技术，不同于常规地震，VSP是在地表附近激发，在井中不同深度布置检波器进行观测。
大地电磁occam反演流程
大地电磁occam反演流程下载温馨提示:该文档是我店铺精心编制而成，希望大家下载以后，能够帮助大家解决实际的问题。
文档下载后可定制随意修改，请根据实际需要进行相应的调整和使用，谢谢!并且，本店铺为大家提供各种各样类型的实用资料，如教育随笔、日记赏析、句子摘抄、古诗大全、经典美文、话题作文、工作总结、词语解析、文案摘录、其他资料等等，如想了解不同资料格式和写法，敬请关注!Download tips: This document is carefully compiled by theeditor. I hope that after you download them,they can help yousolve practical problems. The document can be customized andmodified after downloading,please adjust and use it according toactual needs, thank you!In addition, our shop provides you with various types ofpractical materials,such as educational essays, diaryappreciation,sentence excerpts,ancient poems,classic articles,topic composition,work summary,word parsing,copy excerpts,other materials and so on,want to know different data formats andwriting methods,please pay attention!大地电磁occam反演流程一、OCCAM反演流程概述OCCAM反演是一种广泛应用于大地电磁(MT)数据解释的稳定收敛的反演方法。
构造模型正演技术在下扬子区初步应用效果
构造模型正演技术在下扬子区初步应用效果邓玉胜;苏树桉【摘要】针对圈闭评价中构造形态的可信度这一问题,笔者选用有限元地震波模拟正演技术,对下扬子区的HQ、DY、BP等典型构造进行模型正演计算.结果表明,它有效地克服了解释过程中的片面性,为落实构造、评价圈闭提供了依据.【期刊名称】《物探与化探》【年(卷),期】2004(028)005【总页数】3页(P410-412)【关键词】圈闭评价;构造形态;构造模型;有限元地震波模拟;下扬子区应用效果【作者】邓玉胜;苏树桉【作者单位】中石化集团,新星石油公司华东分公司,规划设计研究院,江苏,南京,210031;中石化集团,新星石油公司华东分公司,规划设计研究院,江苏,南京,210031【正文语种】中文【中图分类】P631.4下扬子区海相中生代、古生代历经多期构造运动，形成的构造圈闭类型多。
圈闭的成因机制分析及评价在很大程度上取决于构造分析，而构造分析可信度又取决于地震测线的控制程度、剖面质量和解释人员的技术水平。
笔者对下扬子区的几个典型构造用有限元地震波模拟系列软件系统进行模型正演计算，从而在提高地震解释水平，验证有关构造圈闭的可信度方面获得了认识。
1 正演程序的原理有限元地震波模拟系列软件系统的核心是用有限元法解二维声波、弹性、粘声和粘弹方程，使用这些解方程程序能够模拟地震勘探中的多种记录剖面(包括单炮记录、时间剖面、偏移剖面、VSP剖面)，将模拟结果与原地震剖面进行对比，拟合程度的好坏影响到构造解释的可靠程度。
构造模型模拟的选用过程就是求解声波方程的过程式中,ρ为密度;α=λ+2μ(λ、μ为拉梅系数);θ=+,其中u=(u,v)位移向量为水平坐标，x为深度坐标,t为时间,f为震源。
2 典型构造圈闭的验证2.1 HQCO2气藏该气藏是下扬子区海相领域的非烃气藏之一，解剖该气藏的圈闭条件及其形成机制对下扬子区海相领域油气勘探具有指导意义。
该气藏海相地层层系发育完整，CO2气田由若干个气藏组成，大多数气藏位于二叠系栖霞组、石炭系船山组和黄龙组、泥盆系五通组海相碳酸盐岩和砂岩地层中。
声波介质一次散射波场高斯束Born正演
声波介质一次散射波场高斯束Born正演岳玉波;钱忠平;张旭东;王德营;岳媛媛;常稳【摘要】Born正演是一种常用的地震波场正演模拟方法,也是线性化地震反演的理论基础.在实际应用时,Born正演通常结合常规的地震射线方法进行实现.为了克服常规地震射线方法的弊端,并且保证地震波场的模拟精度和计算效率,本文提出了一种基于高斯束的一阶散射波场Bo rn正演方法.该方法分为两个环节:首先,我们利用高斯束的走时和振幅信息将地下散射点处的反射率映射为地表束中心位置处的局部平面波;然后,我们利用逆倾斜叠加将局部平面波转化为接收点处的时空域散射波场.在具体的实施过程中,我们提出一种以w avelet-bank方式实现的局部平面波合成方法,同现有的算法相比,可以在保持计算精度的同时,大大减少计算时间;此外,我们还利用最速下降法优化了高斯束的迭代循环过程,进一步提高了Born正演的计算效率.两个模型的应用效果证明,本文所提出的高斯束Born正演方法可以精确、高效的实现声波介质一次散射波场的正演模拟,为三维大规模地震波场的正演问题提供了一种切实可行的实现方案.【期刊名称】《地球物理学报》【年(卷),期】2019(062)002【总页数】9页(P648-656)【关键词】一次散射;高斯束;Born近似;声波介质【作者】岳玉波;钱忠平;张旭东;王德营;岳媛媛;常稳【作者单位】东方地球物理公司物探技术研究中心,河北涿州 072751;山东科技大学山东省沉积成矿作用与沉积矿产重点实验室,山东青岛 266590;东方地球物理公司物探技术研究中心,河北涿州 072751;东方地球物理公司物探技术研究中心,河北涿州 072751;山东科技大学山东省沉积成矿作用与沉积矿产重点实验室,山东青岛266590;东方地球物理公司物探技术研究中心,河北涿州 072751;东方地球物理公司物探技术研究中心,河北涿州 072751【正文语种】中文【中图分类】P6310 引言地震波场的正演模拟是勘探地震学的重要研究课题，其不但可以用来认识地下介质的波场传播规律，进而指导地震资料的采集、处理和解释等环节，还可以利用其共轭算子(或者逆算子)来进行地下介质成像和参数反演运算(Tarantola,1984; Bleistein,1987; Beylkin and Burridege,1990; Chapman,2004).目前，地震勘探规模越来越庞大，所面临的地下地质构造也愈加的复杂，在这种背景下发展具备高精度和高效率的波场正演方法具有重要的理论意义和应用价值.现有的地震波场正演模拟方法主要包括两类：一类是波动方程的数值求解方法，包括有限差分法、有限元法、频谱法等(Tessmer, 2000; Virieux et al., 2011).该类方法直接对波动方程进行数值求解，可以计算直达波、折射波、一次反射波和多次波等全部的波场信息，具有很高的波场模拟精度.但是该类方法的计算效率较低，特别是在模拟高频地震波场或者模型速度较低时，为保证计算过程的稳定性而采用的小计算网格往往会导致计算量大幅度增长.因此，较低的计算效率严重地制约了波动方程数值解法在三维大规模地震波场正演时的应用.另一类是射线类方法1983, 2005; 李辉等, 2012)，该类方法以波动方程的高频近似解为基础，通过射线追踪来逐步求取每个射线节点的走时和振幅，进而进行一次反射地震波场的合成.常规的射线类算法(如两点射线追踪)，计算效率较高，但是其对模型的构造细节不敏感，且存在波场焦散区、阴影区等固有缺陷，致使其波场模拟精度较低.此外，在常规的射线方法中，地下的反射界面是通过离散的样点并结合样条函数进行描述的，当地下构造复杂时，准确的界面描述也是一项很大的难题.射线Born近似或Kirchhoff近似也是常用的地震波场正演方法，具有优于常规射线方法的波场模拟精度(Hudson and Heritage, 1981; Cohen et al., 1986; Jaramillo and Bleistein, 1999; Bleistein et al., 2001)，但是该类方法依赖常规的射线追踪来求取地下走时和振幅信息，因此仍然受到常规射线方法局限性的制约，而且在模拟多次走时波场时存在较大的实现难度.针对上述问题，Huang等(2016)提出了一种基于高斯束的地震波场Born正演方法，该方法使用高斯束取代常规的地震射线进行地下散射点走时和振幅信息的计算，在理论上克服了常规射线方法的弊端.他们在实现过程中使用最小虚值走时点来构建等时线，进而通过等时线叠加求取模拟波场.由于高斯束表征的格林函数是一系列高斯束的叠加，仅仅选择最小虚值走时点显然无法代表全部的波场信息，致使文中给出的高斯束Born正演模拟结果同波动方程有限差分正演结果存在较大的差异.此外，他们所提出的算法是逐道进行的，且在频率域中实现，因此计算效率也不存在优势.本文提出了一种新的声波介质高斯束Born正演模拟算法.该方法将Born正演过程分解为两个基本环节：首先，我们利用高斯束的走时和振幅信息，将地下散射点处的反射率函数映射为地表稀疏束中心位置处的局部平面波；然后，我们通过逆倾斜叠加将局部平面波累加到接收点波场中，得到最终的波场正演结果.在实现过程中，我们使用了Hill(2001)和Gray(2005)所提出的最速下降法来减少循环迭代次数，在保证模拟精度的同时，显著地提高了计算效率.此外，我们在平面波合成时，给出了一种wavelet-bank算法，相对于传统的wave-packet算法计算效率提高了数倍.在下文中，我们首先对相关的理论公式进行推导，接下来对具体的计算实现细节进行探讨，最后通过两个数值模型对波场模拟精度和计算效率进行验证.1 方法原理1.1 基于高斯束的Born正演公式在各向同性声波介质中，假设震源和接收点的空间位置矢量分别为xs和xr，那么由震源xs激发，途径地下散射点x，最后到达接收点xr的一次散射地震波场u(xr,xs,ω)可以利用Born正演公式(Tarantola, 1984; Beydoun and Mendes, 1989)来表示：(1)式中，D为地下散射点的集合，G(x,xs,ω)为由震源下行格林函数，G(x,xr,ω)为接收点上行格林函数，F(ω)为频率域的震源函数，m(x)≈2c1(x)/(x)为地下散射点x 处速度扰动c1(x)所引起的地震反射率，c0(x)为该点的背景速度.为了保证格林函数的正则性，并且使正演过程具备处理多次走时波场的能力，我们使用高斯束来计算式(1)中的格林函数，其表达形式为一系列具有不同出射方向的高斯束的叠加(Hill，2001)：∬×exp[iω τGB(x,x′)],(2)其中，p=(px,py,pz)为中心射线在出射点x′处的初始射线参数，AGB和τGB分别为高斯束的复值振幅和走时，并且具有如下的形式：AGB(x,x′)=AR(x,x′)+iAI(x,x′),τGB(x,x′)=τR(x,x′)+iτI(x,x′),(3)式中，下标R和I分别代表复数参量的实部和虚部.将式(2)、(3)代入式(1)，可以得到以高斯束表征的Born正演公式：∬∬×[AR(xr,xs)+iAI(xr,xs)]exp[iω τR(xr,xs)-ω τI(xr,xs)],(4)其中，AR(xr,xs)和AI(xr,xs)为地下散射点处震源和接收点高斯束振幅乘积的实部和虚部，τR(xr,xs)和τI(xr,xs)为相应的高斯束走时和的实部和虚部，具有如下形式：AR(xr,xs)=AR(x,xr)AR(x,xs)-AI(x,xr)AI(x,xs),AI(xr,xs)=AR(x,xr)AI(x,xs)+AI(x,xr)AR(x,xs),τR(xr,xs)=τR(x,xr)+τR(x,xs),τI(xr,xs)=τI(x,xr)+τI(x,xs).(5)直接利用式(4)进行波场的计算，需要在每个接收点位置计算高斯束.为提高计算效率，我们选择在稀疏的地表束中心位置计算高斯束，并且通过引入一个相移校正量exp[-iωpL·(xr-L)]，来近似计算接收点上行格林函数(Gray and Bleistein, 2009；岳玉波等，2012)∬×exp[iω τGB(x,L)-iωpL·(xr-L)\,(6)并且为了减小上述近似的误差，我们在式(4)中引入下述的单位分解公式(岳玉波，2011)(7)得到最终的Born正演公式：∬dpLxdpLyU(L,xs,pL,ω)(8)其中，为束中心间隔，ωr和w0分别为所选择的参考频率和高斯束初始宽度.U(L,xs,pL,ω)为束中心处L所合成的射线参数为pL的局部平面波，具有如下表达形式：∬×[AR(L,xs)+iAI(L,xs)]×exp[iω τR(L,xs)-ω τI(L,xs)].(9)式(8)的作用是将每个束中心处所合成的局部平面波累加到临近的接收点波场中，该式实际上是一个局部逆倾斜叠加过程，因此可以在频率波数域进行高效计算.式(9)的作用是通过高斯束的振幅和走时信息，将地下的反射率映射为束中心处的局部平面波.该式是本文所提出的Born正演方法的核心，其计算过程决定了最终的波场模拟精度和计算效率，因此我们接下来对其进行重点讨论.1.2 多波至模式直接求解式(9)，需要对震源和束中心高斯束的所有组合对进行循环计算.在三维情况下，若psx,psy,pLx,pLy的采样数分别为Np，那么需要Np4次的束循环运算(我们称该实现方式为全波至模式)，计算量庞大.为了减少计算量，我们应用Hill(2001)，Gray(2005)在高斯束偏移中应用的最速下降法对式(9)的计算进行简化，具体的实现过程如下：(1) 首先，我们通过下式将震源和束中心射线参数转化为中心点和偏移距射线参数pm=ps+pL,ph=pL-ps,(10)得到以中心点和偏移距射线参数索引的高斯束组合；(2) 接下来，对于中心点射线参数为的高斯束组合，按照ph进行循环找到使地下散射点处虚值走时最小的高斯束对.假设此时的偏移距射线参数为则该高斯束对中震源高斯束的射线参数为束中心高斯束的射线参数为计算(11)式中的积分，并且将计算结果累加到局部平面波波场中×exp[iω τR(L,xs)-ω τI(L,xs)],(11)式中，为被震源和束中心高斯束同时照明到的地下散射点，行列式detT为应用最速下降法后产生的振幅校正项(Gray and Bleistein, 2009)detT=(12)其中，Qs,QL分别为震源和束中心高斯束所求取的2×2复值动力学射线追踪参数矩阵(Hill, 2001; Gray and Bleistein, 2009)，c0(x),c0(L)分别为震源和束中心处的速度.(3) 最后，重复步骤(2)直到所有的中心射线参数处理完成，得到对应不同的接收点射线参数的局部平面波U(L,xs,pL,ω).在计算式(8)后，即可得到最终的Born正演结果.上述实现过程可以称为多波至模式，其可以对地下的大部分波场传播路径进行处理(即使存在多次走时波场)，但是所需的束循环次数却由Np4减少为4Np2次.此外，由于高斯束所携带的走时和振幅信息是平滑变化的，求取最小虚部走时的束循环过程可以在稀疏的粗网格上高效的进行.多次波至模式具有比全波至模式高得多的计算效率，而计算精度却损失不大，我们将在模型试算中对此进行验证.1.3 Wavelet-bank算法在将束循环过程进行简化后，我们接下来探讨如何高效的实现式(11)所描述的局部平面波的合成过程.由于式(11)的计算位于整个正演模拟过程的核心位置，其求解过程在很大程度上决定了整个高斯束Born正演的精度和效率.接下来，我们首先对现有的频率域算法和wave-packet算法进行简要的介绍，然后给出一种以wavelet-bank方式实现的时间域高效算法.为了公式表达的简洁，我们将分别用AR,AI来代表式(11)中振幅项的实部和虚部，用τR,τI来代表实部和虚部走时. (1) 频率域算法该算法直接在频率域进行计算，然后利用反傅氏变换求取时间域的局部平面波.由于没有使用任何近似，因此其计算精度是最高的.然而，由于我们需要对每个频率成分进行一次式(11)的循环求解，因此计算量巨大.(2) Wave-packet算法由于频率域算法效率很低，提出了一种wave-packet算法，该方法将频率域的振荡积分转换为了时间域褶积的形式，我们在此直接给出式(11)应用wave-packet 算法后的表达公式(13)其中，是时间域的震源子波，fH(t)为f(t)的希尔伯特变换.由于这两个震源子波项是独立的，我们只需要将地下散射点的反射率通过高斯束的走时和振幅映射叠加到时间域地震道中，待所有的循环计算完成后，将时间域地震道同震源子波项进行褶积即可得到合成的平面波波场.该算法虽然简单直接，但是我们需要将地下单个散射点的反射率m(x)映射到时间域地震道中τR周围的多个时间样点(取决于式中分数项的衰减程度.例如，当τI=5 ms时，时域地震道中|t-τR|≈20 ms的时间样点的振幅项才会衰减为最大值的百分之一).由于散射点的循环位于Born正演的最内层，因此，虽然wave-packet算法的计算效率优于频率域算法，但仍难以令人满意. (3) Wavelet-bank算法我们在此提出一种基于wavelet-bank方式的时间域高效算法.首先，将式(11)通过反傅氏变换转化为时间域：(14)其中，(15)为应用虚部走时对子波频谱进行衰减后的时间域震源子波，为其希尔伯特变换.式(14)只需要将反射率m(x)映射到时间域地震道对应τR的唯一时间样点处，但是每完成一个点的映射我们就需要一次褶积运算，计算效率很低.然而，我们发现式(15)中与虚部走时相关的振幅衰减项只会影响子波的振幅，不会改变子波的相位，因此，我们可以根据所允许的最大虚部走时=5/ωr，创建一个规则采样的虚部走时序列：τI(n)=nΔτI, 0≤n<N,(16)其中，N为离散样点数，ΔτI=/N为采样间隔.然后计算出对应每个序列点τI(n)的震源子波：(17)这样，我们便可以利用线性插值将式(15)表示为：(18)其中，τI(n)=int(τI/ΔτI)×ΔτI，int为取整函数.此时，式(14)便可以表示为：(19)其中， Rn(t)，In(t)为：(20)为满足如下条件的地下散射点的集合：(21)W为线性加权系数，根据虚部走时的大小，其取值如下：(22)Wavelet-bank算法的实现过程可以大致概括为：首先，对于地下每一个散射点，根据其对应的实部和虚部走时计算式(20)，从而将反射率加权映射到对应的Rn(t)，In(t)中；然后，在所有散射点计算完成后，计算式(19)即可得到最终模拟到的局部平面波u(L,xs,pL,t).随着虚部走时的增大，式(15)中的震源子波是单调且平滑衰减的，因此我们不需要选择太大的虚部走时离散样点数N.在实际计算中，我们发现N=10即可完全满足精度的要求.应用wavelet-bank算法，我们只需对地下每个散射点进行两次数据累加运算，此后的褶积运算可以在频率域内高效进行.因此，wavelet-bank算法的计算效率要远高于频率域算法和wave-packet算法.在接下来的模型试算中，我们同样对此进行测试和验证.2 模型试算在本节中，我们将通过两个模型数据进行高斯束Born正演的测试.首先使用一个层状介质模型来对比wavelet-bank算法和现有方法的效率和精度.接下来，使用Marmousi模型测试高斯束Born正演在复杂模型中的应用效果.我们分别计算对应全波至和多波至模式的波场正演结果，并且同具备时间2阶、空间8阶精度的交错网格波动方程有限差分正演结果进行对比.2.1 层状介质模型首先，使用一个层状模型测试wavelet-bank算法的计算效率和精度.图1a展示了所使用的速度模型，该模型横向采样点数为1000，采样间隔为10 m，纵向采样点数为550，采样间隔为5m.图1b为对速度场进行平滑后计算得到的反射率剖面.我们将炮点放在模型表面水平位置5.0 km处并使用主频为20 Hz的雷克子波作为震源，240个接收点均匀的分布在炮点两侧，间隔为20 m(如图1a所示)，每个接收道道长为2 s，时间采样间隔为2 ms.我们利用多波至模式进行高斯束Born正演计算，并且在合成局部平面波时分别使用了频率域算法、wave-packet算法和wavelet-bank算法.我们将频率域计算结果作为基准，来比较另外两种算法的精度和效率.图2a、2b、2c分别为应用频率域算法、wave-packet算法和wavelet-bank算法计算得到的单炮正演记录，可以看到三者几乎没有任何的区别.为了进一步对比精度，我们将频率域算法结果d0作为精确解，然后使用函数ε=‖d-d0‖2/‖d0‖2来计算wave-packet或wavelet-bank算法结果d相对于d0的误差.在对比精度的同时，我们还比较了三种方法的计算时间，最终的时间和精度对比结果如表1所示，可以看到wavelet-bank和wave-packet算法的计算精度类似，同频率域算法的相对误差均在1%以内，但是wavelet-bank算法的计算效率却比其他两种方法高得多，大约是wave-packet算法的7倍，是频率域算法的100倍以上.因此，wavelet-bank算法在保证计算精度的同时，极大提高了计算效率.图1 层状介质模型(a) 速度场，地表黑色线段标识了接收点位置，炮点位于接收道中心； (b) 反射率剖面.Fig.1 Multi-layer model(a) The velocity field, wherethe black line segment on the surface marks the recording receivers andthe shot is at the center of the receivers; (b) The corresponding reflectivity profile.图2 应用不同算法得到的单炮记录(a) 频率域算法； (b) Wave-packet算法； (c) Wavelet-bank算法.Fig.2 Shot records simulated with different modeling algorithms(a) Frequency approach; (b) Wave-packet approach; (c) Wavelet-bank approach.在进行波场模拟时，我们采取了高斯束偏移所使用的束宽度、射线参数采样间隔等参数选取准则(Hill, 2001; 岳玉波, 2001).同高斯束偏移类似，高斯束Born正演方法的波场模拟效果对于上述参数的选取并不敏感，具有较高的冗余度.表1 不同算法的精度和效率对比Table 1 The comparison of accuracy and efficiency between different modeling algorithms频率域算法Wave-packet算法Wavelet-bank算法计算时间(s)127.347.331.05相对误差(%)0.000.970.89 2.2 Marmousi模型接下来，使用Marmousi模型验证高斯束Born正演在复杂模型中的应用效果.图3展示了Marmousi模型的速度场，模型横向采样点数为369，采样间隔为20 m，纵向采样点数为375，采样间隔为8 m，相应的反射率剖面如图4所示.在正演时我们同样使用20 Hz的雷克子波作为震源函数，将炮点放置在横向位置为4.4 km 处的地表位置(如图3所示).共有201个接收道，道间隔为20 m，且沿炮点对称分布.可以看到，该炮对地下照明的区域恰好为模型中构造最复杂的部分.我们分别使用全波至和多波至模式并结合wavelet-bank算法进行单炮记录的计算，所得到的结果如图5a和5b所示，可以看到两者的模拟效果基本相同.为了更好的验证高斯束Born正演的波场模拟精度，我们使用有限差分法进行了相同位置处单炮记录的计算，所得到的结果如图5c所示.可以看到，高斯束Born正演结果同有限差分正演结果的波场整体分布特征相似，证明了本文方法在进行复杂模型正演时的有效性.然而，由于Born正演只能模拟一次散射波场，而有限差分法模拟了地下全部的波场信息，因此两者在局部细节上存在一定的差异.对三者的计算时间进行统计，全波至模式为3.89 s、多波至模式为0.67 s、有限差分法为31.29 s，多波至模式具有明显的效率优势.为对比全波至和多波至的差异，我们计算了两者的残差(图5d)，并求得两者的相对误差为9.7%.相对于全波至模式，多波至模式结果精度略有损失，但是两者的差异主要在于最速下降法近似对波场整体背景振幅的改变，同时考虑到多波至模式的具有接近于全波至模式6倍的计算效率，因此我们认为上述精度损失完全可以接受.为了进一步验证高斯束Born正演对于复杂模型的模拟效果，我们应用多波至模式模拟了总共135炮的地震记录，并且使用高斯束偏移对模拟记录进行了成像处理，所得到的偏移剖面如图6所示.可以看到，最终的偏移结果准确的恢复了该模型复杂的构造形态，进一步验证了本文所提出的高斯束Born正演算法的正确性和有效性.图3 Marmousi模型速度场地表黑色线段标识了接收点位置，炮点位于接收道中心.Fig.3 The velocity field of Marmousi modelWhere the black line segment on the surface marks the receivers and the shot is at the center of the receivers.图4 Marmousi模型反射率剖面Fig.4 The reflectivity profile of Marmousi model图5 应用不同正演方法得到的单炮记录(a) 全波至模式高斯束Born正演； (b) 多波至模式高斯束Born正演； (c) 有限差分法； (d) 全波至模式(图a)同多波至模式(图b)结果之差.Fig.5 The simulated shot gathers using different modeling methods(a) Full arrival Gaussian beam Born modeling; (b) Most arrival Gaussian beam Born modeling; (c) Finite-difference modeling; (d) The difference between (a) and (b).图6 多波至模式模拟数据的高斯束偏移结果Fig.6 Gaussian beam migration image produced by the simulated data with most arrival mode3 结论本文提出了一种兼具模拟精度和计算效率的声波介质高斯束Born正演方法.该方法使用高斯束作为波场传播算子，不但有效的克服了常规射线方法所具有的阴影区、焦散区等固有缺陷，还具备了模拟多次走时波场的能力，可以实现复杂介质一阶散射地震波场的精确模拟.与此同时，本文还着重探讨了该方法的具体实现过程，给出了一套切实可行的实现方案，在保证波场模拟精度的前提下，极大的提高了计算效率.接下来，我们将进一步探讨弹性各向同性介质中地震波场的Born正演方法.此外，我们还将以本文算法为基础构建相应的共轭算子，进而开展最小二乘高斯束偏移技术的相关研究.References【相关文献】Beydoun W B, Mendes M. 1989. Elastic ray-Born L2-migration/inversion. Geophysical Journal International, 97(1): 151-160.Beylkin G, Burridege R. 1990. Linearized inverse scattering problems in acoustics and elasticity. Wave Motion, 12(1): 15-52.Bleistein N. 1987. On the imaging of reflectors in the earth. Geophysics, 52(7): 931-942. Bleistein N, Cohen J K, Stockwell J W Jr. 2001. Mathematics of Multidimensional Seismic Imaging, Migration, and Inversion. New York: Springer.V. 1983. Synthetic body wave seismograms for laterally varying layered structures by the Gaussian beam method. Geophysical Journal International, 73(2): 389-426.V. 2005. Seismic Ray Theory. Cambridge, UK: Cambridge University Press.Chapman C. 2004. Fundamentals of Seismic Wave Propagation. Cambridge, UK: Cambridge University Press.Cohen J K, Hagin F G, Bleistein N. 1986. Three-dimensional Born inversion with an arbitrary reference. Geophysics, 51(8): 1552-1558.Gray S H. 2005. Gaussian beam migration of common-shot records. Geophysics, 70(4):S71-S77.Gray S H, Bleistein N. 2009. True-amplitude Gaussian-beam migration. Geophysics, 74(2): S11-S23.Hill N R. 2001. Prestack Gaussian-beam depth migration. Geophysics, 66(4): 1240-1250. Huang X G, Sun H, Sun J G. 2016. Born modeling for heterogeneous media using the Gaussian beam summation based Green′s function. Journal of Applied Geophysics, 131: 191-201.Hudson J A, Heritage J R. 1981. The use of the Born approximation in seismic scattering problems. Geophysical Journal International, 66(1): 221-240.Jaramillo H H, Bleistein N. 1999. The link of Kirchhoff migrationand demigration to Kirchhoff and Born modeling. Geophysics, 64(6): 1793-1805.Li H, Feng B, Wang H Z. 2012. A new wave field modeling method by using Gaussian Packets. Geophysical Prospecting for Petroleum (in Chinese), 51(4): 327-337.Tarantola A. 1984. Inversion of seismic reflection data in the acoustic approximation. Geophysics, 49(8): 1259-1266.Tessmer E. 2000. Seismic finite-difference modeling with spatially varying time steps. Geophysics, 65(4): 1290-1293.Virieux J, Calandra H, Plessix R É. 2011. A review of the spectral, pseudo-spectral, finite-difference and finite-element modelling techniques for geophysical imaging. Geophysical Prospecting, 59(5): 794-813.Yue Y B. 2011. Study on Gaussian beam migration in complex medium [Ph. D. thesis] (in Chinese). Qingdao: China University of Petroleum.Yue Y B, Li Z C, Qian Z P, et al. 2012. Amplitude-preserved Gaussian beam migration under complex topographic conditions. Chinese Journal of Geophysics (in Chinese), 55(4): 1376-1383, doi: 10.6038/j.issn.0001-5733.2012.04.033.附中文参考文献李辉, 冯波, 王华忠. 2012. 波场模拟的高斯波包叠加方法. 石油物探, 51(4): 327-337.岳玉波. 2011. 复杂介质高斯束偏移成像方法研究[博士论文]. 青岛: 中国石油大学.岳玉波, 李振春, 钱忠平等. 2012. 复杂地表条件下保幅高斯束偏移. 地球物理学报, 55(4): 1376-1383, doi: 10.6038/j.issn.0001-5733.2012.04.033.。
基于三维正演的音频大地电磁阻抗相位不变量校正技术
基于三维正演的音频大地电磁阻抗相位不变量校正技术阮帅;张炯;孙远彬;王绪本【摘要】复杂地形、地质条件的大地电磁数据解释容易出现假象,采用三维正演技术模拟地形和地表不均匀体的背景响应,对实测数据阻抗相位不变量进行校正,实现更准确的定性分析;对三维异常体模型的合成数据进行一维、二维多参数反演试算,以确定地形剧变区选择反演技术的最佳方案.合成数据的试反演结果显示一维反演水平切片假异常较多,二维反演能压制测向假异常,但不能压制走向的假异常,水平切片多出现测向条带.使用本文提出的阻抗相位不变量校正法扣除地形、地表背景响应,结合一维、二维反演,能使实际资料解释成果更加可靠.【期刊名称】《地球物理学报》【年(卷),期】2015(058)002【总页数】12页(P685-696)【关键词】音频大地电磁;阻抗相位不变量;交错网格有限差分正演;Occam反演;共轭梯度反演【作者】阮帅;张炯;孙远彬;王绪本【作者单位】油气藏地质及开发工程国家重点实验室(成都理工大学),成都610059;西安交通大学人居环境与建筑工程学院,西安 710054;贵州省地质矿产勘查开发局103地质大队,贵州铜仁554300;油气藏地质及开发工程国家重点实验室(成都理工大学),成都610059【正文语种】中文【中图分类】P6311 引言音频大地电磁法（以下简称AMT）在我国油气、矿产勘探领域广泛应用（何展翔等，2004）.由于使用的设备轻便、勘探成本低廉，相对于其他电磁法更适应在西部山区开展工作，但山区的剧变地形和地表不均匀性会给数据解释带来很多假象（徐利明和聂在平，2005；赵国泽等，2007）.大量的二维数据正演模拟结果表明，视电阻率、阻抗相位及其旋转不变量不容易受静态效应的影响，因此可作为实际资料一、二维反演结果好坏的重要评价手段（陈小斌，2010；阮帅，2013；苗景春等，2013）.随着三维正反演技术近年来的迅猛发展，更复杂的模型的三维模拟变得可行（Mackie等，1993；Zhdanov和Fang，1996；Siripunvaraporn等，2002；鲁来玉等，2003；沈金松，2003；汤井田等，2007，2008），大量正演结果表明三维的强静态效应也会造成阻抗相位的畸变.因此，和二维模式不同，复杂地形、地质条件下的三维阻抗相位及旋转不变量必须进行有效校正才能用于异常定性分析.虽然近年来三维反演问题的研究进步飞速（Newman和 Alumbaugh，2000；Siripunvaraporn和Egbert，2001；胡祖志等，2006；汤井田等，2007），但研究对象仍基本为水平地形模型的合成数据或大尺度的缓地形工区数据.山区高密数据集的大地电磁三维反演问题及其成效评价仍鲜有研究.笔者认为，在目前AMT 勘探领域的主流定量分析手段仍然是一、二维的情况下，如何有效避免低维反演技术的静位移假象仍是非常值得探讨的问题，这一研究课题的结果可用来评价三维反演效果.本文基于相对成熟可靠的交错网格大地电磁三维有限差分正演，提出简单实用的实测资料阻抗相位不变量校正技术，所有研究基于实际的西南某地金伯利岩AMT勘探的真实地形和真实数据.实测数据的一、二维反演参数基于三维正演合成数据试验选取，以获得低维假设下最合理的反演结果，同时试图结合反演和校正阻抗相位不变量以提高地质解释可靠性.该项目是三维复杂地形、地表不均匀的地质条件下寻找三维地质目标体的典型案例.2 工区概述及三维地电模型我国西南山区地质现象复杂、地形起伏，大部分地区在AMT勘探尺度上（1000m左右探深）几乎没有固定的地层模式.本文讨论的测区地表多出露石冷水组（∈2s）白云岩，在缓坡地带有零星第四系松散堆积物，地表电阻率分布非常复杂，地形变化非常剧烈，1km2范围内的最大高差几乎达到200m，施工测线设计只能顺地形进行，且地下构造非常复杂、无固定走向（不利于二维反演），施工布置见图1.图1 西南某工区高程及测线、测点布置图测线号在顶部标注；三角形为测点，右边数字为测点号.Fig.1 Elevation map with survey lines and sites of the project Line No.on top；site present as triangles，numbers right side are site name.该区勘探目标为金伯利岩管道，实际施工测线距200m，测点距50m，由于地形和植被原因L10线无法施工，实测数据静态效应严重，野外数据初步反演结果出现很多零星异常，且同一测线、不同反演参数和反演方法的结果完全不一致，地质解释人员几乎无法根据实测数据和反演结果得到确切、可信的结论.为寻找准确的金伯利低阻异常、避免三维复杂地质条件下由于不满足一、二维反演应用条件而造成的错误解释，并需将反演结果和实测资料进行互相验证.基于近年MT技术的应用经验，阻抗相位数据不容易受静态位移影响，是非常重要客观的定性参数.但地形剧变的AMT正演研究表明，两个极化方向的阻抗相位和其旋转不变量都被地形影响并无法从中有效提取深部电阻率异常.为此，若能求出地形、不均匀体引起的阻抗相位不变量，然后把这一响应在实测阻抗相位不变量中扣除，则校正后的阻抗相位不变量即为仅包含深部地质体信息的数据.3 基于三维正演的阻抗相位不变量校正技术为探讨如何校正阻抗相位不变量，我们可以利用第2节的地形数据和该区钻井资料建立两个三维地形剧变模型，以探讨地形、地表不均匀对深部异常的影响.第一个模型（MB）为地形、地电背景模型，地表电阻率可填充为观测数据的高频视电阻率或浅层电阻率调查结果，深部电阻率可为均匀空间或典型钻井的分层模型；第二个模型（MA）为异常模型，是在MB背景模型的基础上设置一些低阻管道.根据前期钻井资料，工区主要岩性自上而下依次为牛蹄塘组（∈1n）灰岩；明心寺组（∈1m）黑色页岩夹石英砂岩；清虚洞组（∈1q）白云质灰岩、白云岩及同生角砾岩；石冷水组（∈2s）白云岩等.白云岩厚度约300～500m，深部岩层的电阻率逐渐降低，电性分层大致为类Q型.目标体金伯利岩管道因破碎严重兼富含云母电阻率更低.所有钻井均无电阻率测井数据，实验室测量露头和岩心的电阻率差别巨大，很难获取典型岩石电阻率，针对上述情况，为了尽量获取最佳MB背景模型，只能选取测区比较典型、静态位移较小的测点进行TE模式数据一维反演，得到简化的背景电阻率分层（图2）.图2 测点405视电阻率、阻抗相位曲线（a）及其TE模式一维Occam反演结果（b）图a中实线为Occam反演模型的正演响应；图b中阴影部分为猜测的简化电性分层.Fig.2 1－D Occam inversion results of sites 405with its′apparent resistivity and impedance phase curves（a）Apparent resistivity and impedance phase curves（resulting model response curve present as solid line）；（b）1－D Occam inversion model result（shadowed part present as layers simplified background layers）.MB背景模型的地形和地表电阻率必须和实际工区一致，首先根据栅格化的地形数据建立三维网格，中心标高大于实际高程的网格单元均填充为空气.空气以下的第一层网格单元根据浅地表调查的电阻率分布填充（本项目地形复杂，并无地表电阻率分布数据，且高频10000～1000Hz段的数据质量非常差，故地表填充为均一值）.其他网格单元电阻率按一维分析结果填充（图2b）.这样得到的模型和实际工区的宏观电性大致相同，同时又保留了引起静态效应的剧变地形和浅层不均匀体（图3a、图4a）.MA模型可在MB模型的基础上修改获取，在L04、L12两条测线位置增加两个垂直低阻管道，这两个管道向深部逐渐扩大，图3b、图4b分别显示了MA的三维网格在L04剖面和L12剖面上的二维网格，图5为这两个管道在海拔500m和0m的水平切面位置.对两个模型在东南西北四个边界及深部底边界均使用1.5倍的放大因子等比扩展10个网格和10个空气网格.采用三维交错采样网格有限差分法（Mackieet al.，1993）对 MB和 MA两个模型均进行了1000～10Hz频带与实测数据频率表一样的30个频率的正演计算（因为实测数据大部分测点在10000～1000Hz频段质量较差），获取TE和TM模式的三维阻抗张量，并依此计算阻抗相位不变量.MB 和MA模型的计算结果如图6所示.计算结果表明无论深部有无低阻金伯利岩管道，MB和MA模型的阻抗相位不变量形态相似，基本和地形变化一致，L04线异常被地形响应完全掩盖、L12线有微弱显示（图6中黑实线框所示）若不扣除地形响应，很难直接使用实测数据的阻抗相位不变量判断深部是否存在低阻陡立管道.若已知MB模型的响应，实测阻抗相位不变量的校正非常简单易行.若MB的阻抗相位不变量为PB，实测资料的同频阻抗相位不变量为PO，则简单地形校正公式为计算结果的PC为校正后的阻抗相位不变量，它已经扣除了地形、地表不均匀性（此时MB的浅层网格需根据浅层电阻率调查或高频视电阻率填充）的三维大地电磁响应.图3 三维模型的L4线剖面电阻率分布图（a）MB模型；（b）MA模型.Fig.3 Profile L4of two 3Dresistivity mesh grid slide（a）MB；（b）MA.图4 三维模型的L12线剖面电阻率分布图（a）MB模型；（b）MA模型.Fig.4 Profile L12of two 3Dresistivity mesh grid slide（a）MB；（b）MA.图5 MA模型不同海拔深度切片显示（a）海拔500m；（b）海拔0m.Fig.5 Two slice maps on different elevation of model MA（a）Elv.on 500m；（b）Elv.on 0m.必须指出，由于实际地形建模和三维电阻率网格填充误差，仅在背景选择比较合适的时候PC才有一般阻抗相位的特性.当背景电阻率误差较大时，PC等于45°并不表示电阻率均匀，但较大的PC对应电阻率降低的趋势，反之对应电阻率升高.图7为使用公式（1）对1000Hz，10Hz两个频率的阻抗相位不变量进行校正的结果，可以明显看出地形影响已经消除，异常位置准确.若在实测数据中应用此技术需事先编辑数据，因为存在干扰，校正前必须对阻抗相位不变量做第一象限变换、剔除跳点和平滑插值等操作.测区实际资料的阻抗相位不变量校正分析结果见第6节.其他定性参数如二维偏离度、极化椭圆、电性主轴方向等平面图在三维地形剧变区域也会出现和阻抗相位类似情况，深部信息被地形、地表不均匀性影响，同样可以使用本文的思路对MB模型进行地形和地表不均匀性校正.因篇幅关系，这里不做赘述.4 合成数据的校正阻抗相位不变量配合一维反演解释大地电磁的一维反演技术已经非常成熟，但在二、三维性较严重（TE和TM模式曲线差别极大）的情况下难以应用，以Occam方法为例（Constable和Parket，1987），在做反演之前必须分析选取哪个模式的数据更好，并且电阻率曲线往往要进行平移（静校正）以避免出现“条带现象”.石油勘探的研究对象一般是沉积岩缓变地层，所以可以根据浅层电法勘探的结果（如TDEM）以高频为基准将视电阻率平移（何展翔等，2004），但这样做并不适合火山岩地形剧变地区（阮帅，2013），这些地区本来就可能存在大量陡立地层或电性体，它们在反演断面上的形态也应该是垂直条带（比如本文需要寻找的金伯利岩管道），在实际资料反演中无法确定哪些曲线需要平移，哪些不需要平移.因此，本文对所有数据都不做任何的平移式静位移校正，只选择横向分辨率最佳的模式做一维反演，获取尽量多的垂直条带，然后配合校正阻抗相位不变量排除假异常.选取L04线进行一维Occam反演试验，输入数据分别为TE模式视电阻率和阻抗相位、TM模式视电阻率和阻抗相位、视电阻率和阻抗相位旋转不变量，结果分别如图8—10所示，比较图8、9和10，TM模式数据的一维结果对低阻管道更敏感，同时更容易受地形影响，在真实异常右边又出现了一条假条带；TE模式横向更不灵敏；不变量的结果则是二者的折中.图11是MA模型TM模式一维Occam反演水平切片与10Hz校正阻抗相位不变量的套合.可以看出虽然一维反演切片显示出和地形趋势类似的低阻假象（一般在高地形出现），但结合等值线图可以排除这些假象，两者均显示低阻的地方与真实低阻管道位置一致.通过以上的合成数据分析，可以得知，本测区的实际资料比较适合的一维反演解释方案如下：对TM模式数据剔除飞点，然后进行一维Occam反演，在其水平切片中寻找校正阻抗相位不变量较高，同时电阻率又较低的区域.5 合成数据的校正阻抗相位不变量配合二维反演解释避免静态位移假象更好的反演方法是二维带地形反演.以非线性共轭梯度法（以下简称NLCG）为例（Rodi和 Mackie，2001）：其中，m为反演模型，F（m）为模型的二维正演响应，V为观测数据的协方差矩阵，d为观测数据，L为拉普拉斯算子的矩阵形式，λ为正则化因子，控制模型修正沿着光滑方向还是数据拟合方向，很多学者对这一反演方法进行了深入研究和改进，目前认为影响反演结果的主要因素有输入数据模式、反演网格、正则化因子等. 二维反演的输入数据以TE模式、TM模式或TE、TM模式联合会得到不同的结果，模型研究（苗景春等，2013）表明，TE模式的二维反演对“体”的倾向更好，而TM模式对地层的“形状”反应更好，联合TE和TM模式的数据分辨率更好，但损失倾向信息.本文对三种情况进行了大量试算比较，对于垂直陡立低阻脉，联合TE＋TM的结果更准确，因此用来反演金伯利岩更适合（试验过程从略）.本工区测线为正南北方向，野外观测以正北为坐标系方向，输入数据可不做旋转直接反演，但研究表明在地下结构为三维体而非二维构造时输入数据旋转到最佳电性主轴（陈小斌，2010）的二维反演效果可能更佳.笔者对L12线进行不旋转和SWIFT旋转到最佳电性主轴的反演，结果如图12—15所示.图6 两个模型的500Hz阻抗相位不变量平面图（a）模型 MB；（b）模型MA.Fig.6 Impedance phase invariant slice maps on 500Hz frequency of two models（a）MB；（b）MA.图7 MA模型阻抗相位不变量的两个频率简单地形校正结果（a）1000Hz；（b）10Hz.Fig.7 Corrected impedance phase invariant slice maps on two frequencies of MA图8 L04线TM模式1DOccam反演结果Fig.8 1DTM mode Occam inversionof line L04图9 L04线TE模式1DOccam反演结果Fig.9 1DTE mode Occam inversion of line L04图10 L04线不变量1DOccam反演结果Fig.10 1DInvariant data Occam inversion of line L04图11 MA模型海拔0m综合解释切片阴影图：TM模式1DOccam反演0m切片；等值线：10Hz校正阻抗相位不变量，粗线高值.Fig.11 Interpretation slice map on elv.at 0m Grid－map：slice map of 1DOccam inversion；contours：corrected impedance phase invariant，high value in bold line.图12 L12线TE＋TM模式二维NLCG反演结果（不旋转；正则化因子0.3）Fig.12 NLCG TE＋TM data inversion results of line L12（no rotation；lambda as 0.3）对比图12和图14，图13和图15，可见对于本项目所在的这种地形剧变条件，在不同的正则化因子下，L12测线不旋转的二维NLCG反演的效果更好（图12和图13）.在其他地形、地质条件下，若实际资料的测线不沿正南北，那么旋转方式可能为测线方位角（对应本文而言是0°，即不旋转）或最佳电性主轴.具体是否旋转，必须使用预设近似模型的合成数据试算确定.通过大量试算，最终确定该区更适宜的二维NLCG反演的输入数据是：不经旋转、不经平移法静校正的TE＋TM模式数据；反演浅部网格深度取10Ωm趋肤深度的1／2比较合适（篇幅所限，大量试算结果从略）.图13 L12线TE＋TM模式二维NLCG反演结果（不旋转；正则化因子3）Fig.13 NLCG TE＋TM data inversion results of line L12（no rotation；lambda as 3）采取以上反演方案对MA合成数据的所有测线进行二维反演，其不同深度的水平切片如图16所示，相对于一维反演结果，二维NLCG反演能在测线方向压制地形影响.但是由于各条测线反演独立进行，无法对走向方向进行约束，切片图往往沿测线方向出现条带，易错误解释为南北向断层或岩性分界带.套合校正相位不变量10Hz切片等值线图后，深部100m切片的L06线北段低阻条带假象可避免，同时在L04线和L12线两个较大的非管状体切片形状的两个低阻带上，校正相位不变量有效的圈出了低阻金伯利岩的水平边界.图14 L12线TE＋TM模式二维NLCG反演结果（旋转到最佳电性主轴；正则化因子0.3）Fig.14 NLCG TE＋TM data inversion results of line L12（swift rotation；lambda as 0.3）图15 L12线TE＋TM模式二维NLCG反演结果（旋转到最佳电性主轴；正则化因子3）Fig.15 NLCG TE＋TM data inversion results of line L12（swift rotation；lambda as 3）图16 MA模型的综合解释切片图（a）海拔500m切片；（b）拔100m切片；影像图为NLCG反演电阻率切片，等值线为10Hz校正阻抗相位不变量，粗线高值.Fig.16Interpretation slice maps（a）Elv.on 500m；（b）Elv.on 100m；grid－map：slice map of 2DNLCG inversion resistivity，contours：corrected impedance phase invariant，high value in bold line.6 实测数据的综合解释综合第3、4、5节的研究，实测数据应首先进行剔除飞点处理，然后使用D＋反演技术（Parker，1980）对TE、TM模式的视电阻率和阻抗相位进行圆滑插值.将圆滑插值后的阻抗相位不变量使用公式（1）扣除MB模型合成数据的同频响应（图17）.图17的黑框标识处表明经相位不变量地形校正后，由高地形引起的高相位异常被有效去除，同时西南角的高相位圈闭面积也相应缩小，东南角的地层更趋向于稳定.疑似管状体的位置在北东向对角线上出现3个.为提高解释分辨率，套合电阻率反演结果切片时我们只绘制校正阻抗相位的高值（55°～60°）等值线.一维反演采用第4节模式，二维反演采用第5节模式，因为这两种模式都是依据MA模型反复试算得到的较合理的反演模式，相对于其他方式，它们的结果可靠性应该更强.所有反演输入数据都根据D＋圆滑进行了飞点剔除.一维、二维反演水平切片和校正阻抗相位同频切片等值线的套合结果如图18和图19a所示.图17 实测数据相位不变量500Hz水平切片（a）原始数据；（b）校正后数据；大于60°的奇异值已经剔除.Fig.17 Measured impedance phase invariant slice map on frequency 500Hz（a）Before correction；（b）After correction；Notice：values bigger than 60are removed.图18 实测数据一维OCCAM反演综合解释切片（a）海拔500m；（b）海拔100m；影像图为电阻率切片，等值线为图17b，只绘制55～60°.Fig.18 Interpretation slice maps of measured data 1DOccam inversion（a）Elv.on 500m；（b）Elv.on 100m.grid－map：slice map of inversion resistivity；contours：equals Fig.17right only plot values from 55～60degrees.图19 实测数据二维NLCG反演综合解释图（a）切片海拔300m，做图方式和图18相同；（b）L08剖面，不旋转，正则化因子3.Fig.19 Interpretation slice and profile maps of measured data 2DNLCG inversion（a）Elv.on 300m，same drawing scheme as Fig.18；（b）Profile L08，no rotation，lambda as 3.图18的阻抗相位不变量高值区域显示低阻体主要位于东南部，主对角线上连续出现3个圈闭状高阻抗相位，这些圈闭在西北角并未出现，所以一维反演的西北部低阻异常是目标金伯利岩管道的可能性较小，这和图19a的二维反演电阻率深部切片揭示的地层关系一致.利用同时出现低电阻率、高校正阻抗相位的组合模式判断低阻体更可靠.不仅如此，由于二维反演无东西向约束，在测区中心的高相位圈闭上（图19a）并没有显示出明显的低阻（相对于西南和东北角），但后期的钻井揭示该阻抗相位圈闭确实是深部金伯利岩管道的反映，而这一点在L08反演的剖面图上也有显示.因此校正阻抗相位不变量应该引起重视，即避免大量假异常，也有可能避免因为反演技术局限而丢失异常的情况.依照这种先做三维背景模型MB和异常模型MA，以MB的三维正演得到的阻抗相位不变量为校正量进行校正，同时对MA的正演响应进行一维、二维反演的参数试验，确定最佳反演方案，反演电阻率水平深度切片和校正阻抗相位频率切片配合解释，寻找反演结果显示低阻、同时校正相位不变量显示高值圈闭的区域，获得比较可靠的解释结果，后期钻探验证已经得到证实.7 结论复杂地形、地质条件下，MT／AMT的实测阻抗相位不变量会受地形、地表不均匀性的影响，若根据实际地形和浅地表电阻率建立合适的背景模型（MB），基于三维正演进行阻抗相位校正能够扣除这些影响.在MB的基础上可建立和工区地质体模式一致的假想模型MA，对MA的三维正演合成数据进行一、二维的反演试验，选择比较合理的反演方案，这一方案大部分情况下比较适合实测数据的反演.不同反演技术得到的水平切片有不同的特点，一般一维反演切片有较多零星分布的假异常；二维反演能有效压制同测线上的假异常，但由于各测线做二维反演时互相无约束，在水平切片图上易产生测线方向的条带，只有三维反演才能有效压制各组数据独立反演引起的“数值差”，避免假圈闭和条带.因此，所有一维、二维反演电阻率的水平切片图都应该配合校正阻抗相位不变量等值线进行解释，寻找二者的共同规律，在校正阻抗相位不变量和反演结果矛盾的地方，需要更细致的判断和分析（如图19a），避免丢掉有用的异常.本文研究的解释方法可以推广到西部山前带油气勘探的低频MT数据解释.在三维反演技术普及性应用之前，这种方法可用来提高一、二维反演断面的地质解释准确性，在三维反演广泛应用之后，校正阻抗相位不变量也可用来进行三维反演结果的进一步评价.致谢感谢国家重大科学仪器设备开发专项基金对研究的大力支持，并特别感谢贵州省地质矿产勘查开发局103地质大队为本文的研究提供宝贵的实地数据和相关地质背景介绍.ReferencesChen X B.2010.Effects of topography to magnetotelluric measured data and two dimensional inversion correction with topography（in Chinese）.／／China Geophysics 2010－26th Annual of CGS＆13th Academic Congress of SSC，325.Constable C S，Parker R L.1987.Occam’s inversion：apractical algorithm for generating smooth models from electromagnetic soundingdata.Geophysics，52（3）：289－300.He Z X，Jia J D，Gou L.2004.The role of non－seismic techniques in exploration and development of oil and gas.Petroleum Exploration and Development（in Chinese），28（4）：70－72.Hu Z Z，Hu X Y，He Z X.2006.Pseudo－three－dimensional magnetotelluric inversion using nonlinear conjugate gradients.Chinese J.Geophys.（in Chinese），49（4）：1126－1234.Lu L Y，Zhang B X，Bao G X.2003.Modeling of three－dimensional magnetotelluric response for a linear earth.Chinese J.Geophys.（in Chinese），46（4）：568－575.Mackie R L，Madden T R，Wannamaker P E.1993.Threedimensionalmagnetotelluric modeling using difference equations：Theory and comparisons to integral equation solutions.Geophysics，58（2）：215－226.Miao J C，Ruan S，Zhang Y.2013.The application of the audio magnetotelluric sounding method to the precise interpretation of formal and reverse faults.Geophysical ＆ Geochemical Exploration （in Chinese），37（4）：681－686.Newman G A，Alumbaugh D L.2000.Three－dimensional magnetotelluric inversion using non－linear conjugate gradients.Geophys.J.Int.，140（2）：410－424.Parker R L.1980.The inverse problem of electromagnetic induction：Existence and construction of solutions based on incompletedata.J.Geophys.Res.，85（B8）：4421－4428.Rodi W L，Mackie R L.2001.Nonlinear conjugate gradients algorithm for 2－D magnetotelluric inversion.Geophysics，66（1）：174－187.Ruan S.2013.Discussion of static shift effect in practical AMT data inversion and interpretation（in Chinese）.／／China Geophysics 2013－29th annual of CGS.volume 24：631－632.Shen J S.2003.Modeling of 3－D electromagnetic response in frequency domain by using staggered－grid finite difference method.ChineseJ.Geophs.（in Chinese），46（2）：281－288.Siripunvaraporn W，Egbert G.2001.An Efficient data－subspace inversion method for 2－D magnetotelluric data.Geophysics，65（3）：791－803. Siripunvaraporn W，Egbert G，Lenbury Y.2002.Numerical accuracy ofmagnetotelluric modeling：a comparison of finite difference approximations.EarthPlanetsSpace，54（6）：721－725.Tang J T，Ren Z Y，Hua X R.2007.The forward modelling and inversion in geophysical electromagnetic field.Progress in Geophysics（in Chinese），22（4）：1181－1194，doi：10.3969／j.issn.1004－2903.2007.04.025.Tang J T，Luo W B，Liu C S.2008.Impulse response of seafloor hydrocarbon reservoir model.Chinese J.Geophys. （in Chinese），51（6）：1929－1935，doi：10.3321／j.issn：0001－5733.2008.06.036.Xu L M，Nie Z P.2005.A fast forward algorithm for modelling vector electromagnetic scattering from buried dielectric objects.ChineseJ.Geophys.（in Chinese），48（1）：209－215，doi：10.3321／j.issn：0001－5733.2005.01.028.Zhao G Z，Chen X B，Tang J.2007.Advanced geo－electromagnetic methods in China.Progress in Geophysics （in Chinese），22（4）：1171－1180.Zhdanov M S，Fang S.1996.Quasi－linear approximation in 3－D electromagnetic modeling.Geophysics，61（3）：646－665.附中文参考文献陈小斌.2010.地形对大地电磁观测资料的影响及二维带地形反演校正.／／中国地球物理2010——中国地球物理学会第二十六届年会、中国地震学会第十三次学术大会论文集，325.何展翔，贾进斗，苟量.2004.非地震技术在油气勘探开发中的作用.石油勘探与开发，28（4）：70－72.胡祖志，胡祥云，何展翔.2006.大地电磁非线性共轭梯度拟三维反演.地球物理学。
盾构管片壁后注浆密实度地质雷达正演模拟与实例分析
这些间隙，就要在盾构机推进过程中，保持一定压力不 间断地从盾尾直接向壁后注浆［2］。壁后注浆分为同步注 浆、及时注浆和二次注浆。可见，管片背后注浆的主要目 的有防止地层变形，减少隧道的沉降量，使管片得到部分 稳定，防止管片偏移。并且浆液凝结后具备一定的强度， 提高了隧道的抗渗能力，增加衬砌接缝的防水功能，可以 说管片背后注浆起到了多方面的作用［3］。
- 44 -
Q 质量检测 uality Test
工程质量
第39 卷
27 层和 30 层。
5 结语
1）本文从原结构和加固工程两个方面介绍了已加 固剪力墙结构的检测内容，从抗震构造措施、加固构造 措施、承载力验算三方面列举出鉴定结果。根据鉴定结 果，建议委托方进行楼板全数检测鉴定，有必要时进行 加固修复，为房屋后续正常使用提供技术支持。
本文 通 过 采用时 域 有限差 分（F D T D）方法的 GPRMax3.0 程序建立不同的管片壁后注浆缺陷 3D 模 型，然后分析其正演模拟反射特征，再对比分析地质雷 达在实际盾构隧道管片检测结果，从而提高技术人员对 雷达图像解释的准确性，有助于更好地指导地质雷达在 盾构管片壁后注浆检测的应用［6］。
GPRMax 是爱丁堡大学的 Antonis Giannopoulos 于1996 年推出来的一种基于时域有限差分（FDTD）算 法求解麦克斯韦方程组和理想匹配层（PML）边界吸收 条件的探地雷达正演数值模拟程序［8］，既可以用于对探 地雷达（GPR）进行建模，也可用于对其他电磁波传播进 行建模。
该程序基于 Yee 网格块［9］（见图 1）对空间和时间 连续体进行离散。将麦克斯韦旋度方程转化为差分方 程，通过解差分方程得到微分方程解的近似值［10］。在时 间上迭代求解，该方法还具有较好的稳定性和收敛性， 现成为最重要电磁场数值方法之一。
三维大地电磁正演及反演方法研究现状
三维大地电磁正演及反演方法研究现状摘要：近年来，随着计算机技术和三维电磁模拟技术的发展。
基于积分方程法(IEM)、有限差分法(FDM)和有限单元法(FEM)的三大方法的三维大地电磁正演模拟技术得到了极大的发展。
基于最优化理论的三维大地电磁反演研究也得到了快速发展。
关键词：电磁正演模拟；数值模拟技术；大地电磁反演1 三维大地电磁正演方法研究现状积分方程法(IEM)、有限差分法(FDM)和有限单元法(FEM)是数值模拟技术中的三大方法。
近年来，基于上述方法的三维大地电磁正演模拟技术得到了极大的发展。
在积分方程法中，麦克斯韦方程组被转换为 Fredholm 积分方程，并以此实现对电磁场散射方程的离散，从而得到与待求电场有关的复线性方程组。
该线性方程组的系数矩阵为致密的复数矩阵。
在简单模型的模拟计算中，该方法仅对异常区进行离散，由此得到规模较小的致密系数矩阵，这有利于线性方程组的快速求解。
基于积分方程法在内存消耗、计算速度等方面的优势，该方法在电磁模拟的研究中受到了研究人员的重视。
然而必须指出的是，在复杂地球物理模型中，必须考虑全区域离散化，此时基于积分方程法得到的系数矩阵表现为大规模的致密矩阵，不利于方程组求解。
因此，考虑到对复杂模型模拟计算的适应性问题，认为基于积分方程法的三维 MT 正演技术在反演中的应用具有一定的局限性。
有限差分法发展最为成熟数值计算方法之一，该方法基于差分原理，以节点的差商近似为相应的偏导数，从而得到节点上关于物理场的相关线性方程组。
在电磁场模拟计算中，该线性方程组的系数矩阵为大型稀疏复数矩阵，基于合适的存储和求解方案，可以较快速的对其进行求解。
早在上世纪 60 年代，有限差分法就被用于地球物理场的模拟计算。
进入上世纪90 年代以后，随着交错网格有限差分理论的提出，该方法在地球电磁场模拟研究领域中得到了更为广泛的关注和重视。
交错网格有限差分法在处理内部电磁差异引起的电场与磁场不连续现象等方面具有相当优势，且易于适合编程实现，因而在三维大地电磁场的正演模拟中得到了广泛应用。
分频属性反演与模型正演技术识别火山岩气藏夹层
分频属性反演与模型正演技术识别火山岩气藏夹层李道清;王彬;仇鹏;苏航;杜庆祥;邹明德【摘要】目前火山岩储层的夹层分布特征研究处于探索和初级定性描述阶段,几乎没有定量的区分技术手段,因而,在低渗透薄储层火山岩气藏有效开发过程中,夹层的识别成为一个新的热点和难点.在分析了克拉美丽气田滴西17井区火山岩气藏储层及其夹层特征的基础上,利用分频属性反演与模型正演相结合的方法能够有效识别出火山岩储层中夹层及其展布特征,而且大大提高了地震的纵横向分辨率.将该方法应用于火山岩储层的开发实践中,取得了很好的应用效果,且对于火山岩气藏夹层识别方法的推广及指导油气藏开发具有重要的意义.【期刊名称】《石油天然气学报》【年(卷),期】2014(036)002【总页数】4页(P66-69)【关键词】低渗火山岩;分频属性反演;模型正演;夹层识别;水平井开发【作者】李道清;王彬;仇鹏;苏航;杜庆祥;邹明德【作者单位】中石油新疆油田分公司勘探开发研究院,新疆克拉玛依834000;中石油新疆油田分公司勘探开发研究院,新疆克拉玛依834000;中石油新疆油田分公司勘探开发研究院,新疆克拉玛依834000;中石油新疆油田分公司勘探开发研究院,新疆克拉玛依834000;中国地质大学(北京)能源学院,北京100083;中石油新疆油田分公司新港作业分公司,新疆克拉玛依834000【正文语种】中文【中图分类】P631.44克拉美丽气田滴西17井区火山岩气藏储层经过风化淋滤改造,使其非均质性大大增强,夹层广泛发育,储层预测困难,也使得水平井实施过程中有效避开夹层难度异常大。
因此,提高夹层的识别和分辨能力是火山岩气藏开发的重点和难点,也是气藏有效开发的关键因素。
目前对于火山岩储层的夹层描述还处于探索阶段或者初级定性描述阶段[1～4],几乎没有定量的区分手段和方法[5],也没有可借鉴的成熟技术方法。
因受地震数据体分辨率和信噪比的制约以及地震解释人员不同的知识建构和经验,地震解释预测的圈闭形态、地层深度和储层特征等与实钻情况存在一定的差异,应用地震属性也往往不能很好地解决储层的横向预测的准确性和多解性问题,所以常规地震方法几乎无法识别和预测夹层的展布特征。