2. 江苏省地震局, 南京 210014
2. Jiangsu Earthquake Administration, Nanjing 210014, China
实际的地电结构一般是三维的,在大地电磁中应该采用三维反演获取地下电性结构.但目前三维反演解释由于存在以下三个方面的问题而没有得到推广使用.首先,迄今为止还没有解决复杂三维正演问题的快速算法,这就使得即使采用并行计算,三维反演计算量仍然非常大;其次,已有三维反演算法存在严重的稳定性问题,为了得到稳定的反演结果,我们通常不得不简化三维模型;最后,考虑面积性三维大地电磁勘探费用问题,一般进行大面积地质调查时不采用三维大地电磁勘探,通常只是布设几条单独剖面进行观测,特别是区域性勘探和深部结构研究.因此,目前对于实际资料的解释主要依赖于二维反演技术.但在对实际资料进行二维解释时,会碰到一些比较棘手的问题,比如:将观测数据旋转到构造走向方向进行二维反演,经常会遇到TE、TM和TE+TM反演结果相差非常大的情况[1],此时选择哪种极化模式数据能最大限度地减小三维效应的影响,得到的反演结果更能反映地下结构的主要特征?这是一个非常重要的问题,因为我们后续的电性结构分析和地质解释都是建立在该反演结果基础上的.已有部分研究者注意到这个问题,开展了相关的研究工作,但关于反演数据极化模式的选择却存在非常大的分歧.Wannamaker等[2~4]通过模拟层状介质中的三维低阻异常体的影响,认为TM模式响应受三维效应的影响较小,因此推荐二维近似反演中采用TM模式数据.但其研究主要集中在对反演结果的定性评价上,对三维和二维模型响应数据的差异、如何认识二维模式的近似程度,如何利用观测的TE模式和TM模式的视电阻率和相位资料等问题并未开展相应的研究.陈小斌等[5]研究了均匀半空间中的三维低阻异常体模型,发现对该模型的二维反演中,采用TM模式更易得到真实的模型信息,但其采用三维模型比较简单,结论是否具有普遍性还有待深入研究.Berdichevsky等[6~8]通过曲线形态的分析,认为TM模式响应受导电三维结构的影响较小,而TE模式受高阻三维结构影响较小,因此提倡采用TE+TM模式联合反演.但其对响应曲线受三维效应影响的研究主要是针对一维响应曲线而言的,并没有分析三维响应和二维响应之间的差异,因此其得出的结论也存在一定的疑问.Ledo等[9]研究了三维数据的二维反演中应注意的一些问题,其研究也主要是基于反演结果的分析,并没有给出确定的反演中数据极化模式选择问题的答案.还有一些研究者[10]从纯粹的二维理论模型出发,得出的结论是TE模式对探测高导层敏感,而TM模式对探测高阻层有效.
很明显,上述各研究之间存在较大的分歧.迄今为止,这些分歧还没有得到彻底的解决.尤其是在国内,一般倾向于不单独采用某一个极化模式数据进行反演,而倾向于TE+TM模式的联合反演.为了使这个问题研究得更透彻,本文采用二维、三维模型技术和二维反演技术,进行了大量的正反演计算,通过分析沿走向方向不同延伸长度三维模型和二维模型的响应差异,以及比较不同延伸长度三维模型响应的二维反演结果,从模型响应差异和二维反演结果合理性两个方面对这个问题进行了深入的研究.
2 模型计算和分析本文三维正演采用Mackie的交错网格有限差分程序[11],我们对其进行了一定的改进,将空气层设定修改为完全自适应,并编写了三维正演数据输出的接口部分.同时我们还对其中相应剖面的二维模型进行了二维有限元正演[12]计算,以其反演结果作为比较的标准.二维、三维正演响应的结果被导入到大地电磁可视化系统(MT-Pioneer)[13]中进行后处理和二维反演.研究中设计了大量二维、三维理论模型进行分析计算,由于篇幅有限,以下仅给出其中的几种模型的计算结果.
2.1 高导异常体模型设计了一个复杂的三维模型(见图 1).模型为二维阶梯模型中赋存一个三维高导异常体,其顶部埋深为0.5km,y方向延伸长度为4km,z方向的延伸深度为2km,沿x方向延伸长度是变化的,分别为4、6、8、10、20km和40km,以及无限延伸的二维模型.在0.001~100Hz频段范围内设计频点40个,三维正演网格为56×61×51个.分析选取的测线剖面靠近三维异常体的中心位置,剖面展布方向为y方向,走向方向为x方向.
图 2、3是三维高导异常体不同走向延伸长度计算的三维响应与二维响应的数据差异分布图.按照测线的展布情况,此时三维XY模式相当于二维主轴方向的TE模式,YX模式相当于TM模式.由图可见,当三维异常体的走向-倾向比为1:1(异常体走向延伸4km)时,三维响应和二维响应的差异最为严重.对于TE模式,视电阻率相对误差出现了大片相对误差90%以上的区域,相位差异也有大片区域超过了20°.而对于TM模式,视电阻率相对误差未超过90%,且大部分都在30%以内,相位差异都在4°以内,差异的程度明显低于三维异常体的走向-倾向比为5:1(异常体走向延伸20km)时的TE模式.当三维异常体的走向-倾向比为2:1(异常体走向延伸8km)时,TM模式视电阻率差异已经不明显了,绝大部分区域视电阻率相对误差都在25%以内,相位差异基本都在2°以内,已基本可用来进行二维反演了.而对于TE模式而言,仅当三维异常体的走向-倾向比(异常体走向延伸40km)为10:1时,其视电阻率相对误差才降到30%以下,相位差异降到4°以下,才能用来进行二维反演.由此可见,对于高导异常体模型,TE极化模式对模型二维性的要求远远强于TM极化模式.
图 4、5分别是高导模型的三维响应与对应二维响应视电阻率和相位误差分布直方图.图中横坐标为误差分布,纵坐标为误差频次.从图中可以看到,当三维异常体的走向-倾向比为1:1时,二维TE模式数据与三维响应之间存在较大的差异,而二维TM模式和三维响应的差异要小于前者;随着延伸长度比继续增大,TM模式的数据差异迅速向误差中心的零点集中,误差分布呈高峰状;而TE模式一直到三维异常体走向-倾向比为10:1时,视电阻率误差还在一个较大的范围内分布,不如TM极化模式那样呈现为明显的高峰状.相比较而言,相位数据受三维畸变的影响较视电阻率要小得多,尤其是TM极化模式,甚至在三维异常体的走向-倾向比为1:1时,TM模式相位的数据差异分布就呈现为明显的高峰状,表明这时数据误差已基本集中在零点附近了.
图 6显示的是高导模型的三维响应和对应二维响应差异随三维异常体走向延伸长度的变化情况.这里设定了参与统计数据误差门限值,视电阻率和相位误差的门限值分别为5%和1°,图中统计的都是大于门限值的误差.从图中可以看到,当三维异常体的走向-倾向比为1:1时,二维TE模式和三维响应的视电阻率误差为40%,相位误差大约为9°,而此时二维TM模式和三维响应视电阻率误差为20%,相位误差约为2°.TE模式二维响应和三维响应误差随三维异常体走向-倾向比增大,其减小的趋势非常缓慢;当三维异常体走向-倾向比为5:1时,TE模式的视电阻率误差仍然有30%,相位误差为5°.只有当三维异常体的走向-倾向比为10:1时,二维TE模式和三维响应视电阻率误差才减小到20%,相位误差接近2°.
以上主要从二维和三维模型响应差异的角度进行了分析,结果表明当三维高导体走向-倾向比较小、三维性强时,二维TM模式和三维响应的数据差异要远小于TE模式.随着三维高导体走向-倾向比增大,在构造上愈来愈接近于二维模型假设时,TE模式的数据响应差别仍然比较大,而TM模式的数据响应已经非常接近于相应的二维模型响应了,特别是TM模式相位的差别更小.因此从数据对比角度而言,在三维高导模型下,TM极化模式对模型二维性的要求远远低于TE极化模式,因而TM极化模式数据更适合做二维反演解释.然而,如前所述,二维理论模型的研究结果表明,TM模式对高导模型不敏感,因此,数据上的差异是否可能体现的是对异常体分辨能力的差别,而不是维性引起的差别?以下将对三维高导模型不同走向延伸长度的响应数据进行二维反演,以进一步对这个问题进行验证.反演主要采用非线性共轭梯度法(NLCG)[14],数据模式分别采用单独TE、单独TM和TE+TM联合反演三种,并和相对应二维模型的反演结果进行对比.
图 7是三维高导模型的二维反演结果.从图中我们看到,TE模式数据的反演结果非常不理想,虽然对三维高导体有一定反映,但反演得到的电阻率值和实际电阻率值差别非常大.三维高导体走向方向延伸长度较小时,反演得到虚假的下覆结构.而TM模式的反演结果相对而言要好得多,反演得到的下伏二维阶梯模型的形态和电阻率值大小也较好,结构位置比较准确.TE+TM模式联合反演的结果虽然也可以接受,但当三维高导体走向延伸长度较小时,也存在比较明显的虚假构造,表明TE+TM模式联合反演中,数据受到了二维性较差的TE模式数据的污染.通过反演对比我们可以看到,对于三维高导异常体模型,单独采用TM极化模式数据反演更能得到真实的地电结构形态,尽管在理想二维高导情况下,TE+TM联合反演的结果要好于任何一种单独极化模式的反演结果.
以上关于高导异常体的三维模型研究表明,无论是从数据差异性方面,还是从二维反演结果对比的角度,大地电磁二维反演解释中采用TM极化模式都更为合理.对于高阻异常体模型,是否也有同样的结论呢?Berdichevsky等研究认为TE模式受三维高阻结构影响较小[7].但在Berdichevsky等的研究中,以一维响应作为正常曲线来考察三维畸变的强弱,故其得到的结论值得商榷.为了对这个问题进行进一步的探讨,我们将图 1所示的高导异常体模型修改为高阻异常体模型,如图 8所示,其几何结构与图 1所示的模型完全相同,但高导异常体改成了高阻异常体,围岩的电阻率值也有所修改.同样计算了三维高阻异常体的走向长度分别为4、6、8、10、20km和40km,以及无限延伸的二维模型.结果对比分析选取的测线剖面位置、展布方向、走向方向与2.1节高导模型完全相同.也同样,我们先做二维TE、TM极化模式与相应的高阻异常体三维响应的数据差异对比,然后,进行各种情况反演结果的对比分析.
图 9、10分别是不同走向延伸长度的三维高阻异常体模型计算的三维响应与二维响应的数据差异分布.与2.1节相同,三维XY模式相当于二维主轴方向的TE模式,YX模式相当于TM模式.由图可见,当三维高阻异常体的走向-倾向比为1:1时(异常体走向延伸4km),三维响应和二维响应的差异最为严重.对于TE模式,视电阻率相对误差在较大区域里大于70%,但相位差异相比三维高导模型而言要小,只有小片区域误差在8°左右.TM模式视电阻率相对误差都在60%以下,相位差异8°以下.可以看到,随着三维高阻体走向长度增加,TE模式视电阻率相对误差减小的速度要快于TM模式,当三维异常体走向-倾向比约为2.5:1时(异常体走向延伸10km),TE模式视电阻率相对误差降到30%以下,和TM模式视电阻率误差差不多.三维异常体走向-倾向比为1.5:1时(异常体走向延伸6km),TE模式相位误差和TM模式相位误差绝大部分已经降到6°以下,表明两种极化模式相位响应受三维高阻体影响较小.
图 11、12分别是高阻模型的三维响应与对应二维响应视电阻率和相位误差直方分布图,从图中可以看到,当三维异常体走向-倾向比为1:1时,TE模式视电阻率误差较大,而TM模式视电阻率误差要小于TE模式.随着三维高阻体走向延伸长度增加,TE模式视电阻率误差向误差中心零点集中的速度比TM模式要稍快一些.当三维高阻体走向-倾向比为10:1时,TE模式视电阻率误差几乎全部集中到零点附近,TM模式视电阻率误差也基本集中到零点附近.与高导三维模型类似,相位响应受高阻三维畸变的影响要远远小于视电阻率响应,不论是TE模式还是TM模式,甚至在三维高阻体走向-倾向比为1:1时,两者相位误差分布就呈现为明显的高峰.
图 13是高阻模型的三维响应和二维响应的误差随三维体走向延伸长度的变化情况.同样对数据误差设定了参与统计数据误差门限值,取值与2.1节高导三维模型相同.从图中可以看到,当三维异常体的走向-倾向比为1:1时(走向延伸长度4km),二维TE模式和三维响应的视电阻率误差为30%,而此时二维TM模式和三维响应视电阻率误差为20%,两者相位误差都差不多,约为4°.但随着三维高阻体走向程度的增大,TE模式视电阻率误差迅速减小,在走向-倾向比在2.5:1和5:1之间(约为3:1)时,已和TM模式视电阻率误差基本相同,而当走向延伸长度为40km时,TE模式视电阻率误差甚至已经减小到接近于零,而TM模式视电阻率误差则在10%左右.两种极化模式相位数据误差本身就比较小,基本上都在4°左右变化,表明相位响应受三维高阻体畸变影响较小.
以上有关高阻异常体数据对比分析的研究表明,当三维高阻体的走向-倾向比较小时,二维TM模式和三维响应的差异要小于TE模式,但随着三维高阻体走向-倾向比增大,TE模式与三维模型响应之间的差异迅速减小.当三维高阻体的走向-倾向比约为3:1时,二者与三维响应的差异基本一致.当走向-倾向比继续增加,TE极化模式与三维响应的差异要比TM模式小.值得注意的是,在高阻异常体模型下,无论哪种极化模式,二维响应与三维响应的差异主要体现在视电阻率上,相位的差别一直较小.从数据对比结果来看,对于高阻异常体模型,TE、TM模式响应与三维响应的差异相差不大.理论上看,TE模式可能稍微好一点,因为走向-倾向比只要大于3:1,TE模式的差异性就要比TM模式的好.但是,考虑到野外实际构造体的延伸情况,如果走向-倾向比大于3:1的高阻体异常体很少见,则TM极化模式要好于TE模式.故从数据对比的角度看,两种极化模式各有优势,难分伯仲.
在进行上述数据对比工作的时候,我们忽视了一个问题.在二维TE模式下,高阻体的响应本身就很弱,亦即TE模式对二维高阻体不敏感.因此,这里我们同样需要验证:TE模式与三维响应的差异小究竟是维性因素引起的,还是因为二维状况下,TE模式对二维高阻体不敏感造成的,对此我们进行了各种情况下二维反演计算.反演方法与2.1节相同,反演的数据模式分别采用单独TE、单独TM和TE+TM联合反演三种,同时和对应的二维模型的反演结果进行对比,反演结果如图 14所示.从图中可见,随着异常体走向延伸长度的增大,所有反演结果越来越趋于二维理论模型的反演结果.很明显,TM模式的反演效果最好,从一开始(走向-倾向比1:1)就与二维理论模型的反演结果很接近,TE+TM模式次之,TE模式最差.因此,高阻异常体模型的反演结果再次证明,对于三维模型,TM模式的反演结果仍然明显好于TE模式的结果,在大多数情况下也优于TE+TM联合反演的结果.
通过三维高导模型序列和三维高阻模型序列的计算,我们系统开展了三维模型和相应二维模型在响应数据差异、反演结果两方面对比研究,结合我们所作的其他模型的结果,发现无论是三维高导模型还是三维高阻模型,TM模式反演结果对真实结构的反应都较好,不会引入太多的虚假构造,一般情况下明显好于TE模式的反演结果,大多数情况下还优于TE+TM模式联合反演的结果.
上述模型计算没有考虑数据误差对反演结果的影响,李墩柱等[15]详细研究过大地电磁中数据误差对反演结果的影响,但他们只讨论了TE+TM联合反演结果.我们选取三维高阻模型序列中走向长度为4km的三维模型和走向无限延伸的二维模型,对模拟数据加入5%满足均匀分布的随机误差,误差加入方式与李墩柱等[15]相同.同样对单独TE、单独TM和TE+TM联合反演结果进行了对比,得到的结果也支持我们上述结论.
上述结论与我们原来基于二维理论模型所得到的认识相差很大.二维理论模型研究的结果认为TE模式对高导体敏感,而TM模式对高阻体敏感,因此,在实践中大多倡导TE+TM模式的联合反演,以利用TE、TM极化模式各自的优势,并增加数据约束量.然而,三维模型研究却表明,在三维高导体模型下,TE模式数据受到了三维高导体散射的影响,畸变严重,恰恰不能用于高导体的探测,反而在高阻体探测的时候起的作用可能更大.这个问题只要稍微深入的思考一下,也是不难解释的.引起大地电磁数据畸变的主要是电场.在TE极化模式下,电场沿着构造无限延伸的走向方向,相对于沿剖面变化的横向不均匀体,是连续变化的切向分量,因此不会引起多少畸变.在TM极化模式下,电场穿越沿剖面变化的横向不均匀体,此时,电场将产生突变,并在电性结构的分界面上形成积累的面电荷,从而引起较大的畸变.但是,我们不能据此得到TM模式数据较TE模式数据更差的结论.事实上,因为理论模型是对实际模型的逼近,TM模式的这种畸变与实际地电模型的要求是相一致的,亦即利用TM模式可以模拟实际地电结构在横向上的各种变化,而TE模式所要求的电场不穿越界面这个条件在实际地电模型面前显得过于苛刻.实际情况下,沿构造走向方向也是存在不少的电性不均匀体的,因此电场也是存在较大的变化的,而这种变化,用二维TE极化模式很难模拟.这种情况与观测实际资料的特点是相符合的.在实测数据中,两个方向的视电阻率都会存在较大的跳跃性,表明视电阻率沿剖面连续变化的TE模式基本上不存在.而有关实际资料的二维反演实践也证明,TM模式与TE模式相比,更易于反演拟合,且冗余构造少.
上述反演结果的对比研究都是基于NLCG反演方法,为了对上述结论做进一步的印证,消除反演方法本身对本文研究结论的影响,我们还采用Occam方法[16]对数据进行了反演,结果也表明单独采用TM模式的反演结果更能反映实际地电结构形态,因此不同反演方法本身对本文结论的影响较小.
以上研究还有一个重要的结论是,相位受三维效应的影响较视电阻率小.这不仅使得我们在反演中应适当提高相位的权值,同时还对我们在功率谱的选择编辑、形成迭加功率谱时提出了新的要求.因为我们以往在进行这项工作时,往往重视视电阻率的编辑,以视电阻率曲线的质量为依据进行功率谱的取舍,而忽视了相位数据的重要性.现在需要把这个标准倒过来,应该更重视相位数据的编辑与取舍.
我们的研究结果已经证明,采用TM模式进行二维反演可以较好地重建三维模型信息,而较少地带入冗余构造.由于实际的大地电磁观测数据可以看作是具有一定走向延伸的三维体的响应,因此,基于TM极化模式数据的二维反演可以获得较为合理的电性结构模型.依据上述研究成果,对于实测数据,应优先考虑采用TM模式数据进行二维反演,其次是TM+TE模式,一般不要单独采用TE模式.
致谢中国地震局地质研究所赵国泽研究员在本项研究上给予了多处指点,中国地质大学王家映教授对高阻体部分的研究提出过建设性意见,在此表示衷心感谢.感谢审稿专家提出的宝贵意见和建议.
[1] | 陈小斌, 赵国泽, 马霄. 关于MT二维反演中数据旋转方向的选择问题初探. 石油地球物理勘探 , 2008, 43(1): 113–118. Chen X B, Zhao G Z, Ma X. Preliminary discussion on selecting rotation direction in 2-D MT inversion. Oil Geophysical Prospecting (in Chinese) , 2008, 43(1): 113-118. |
[2] | Wannamaker P E, Hohmann G W, Ward S H. Magnetotelluric responses of three-dimensional bodies in layered earths. Geophysics , 1984, 49(9): 1517-1533. DOI:10.1190/1.1441777 |
[3] | Wannamaker P E, Booker J R, Jones A G, et al. Resistivity cross section through the Juan de Fuca subduction system and its tectonic implications. Journal of Geophysical Research , 1989, 94(B10): 14127-14144. DOI:10.1029/JB094iB10p14127 |
[4] | Banks R J, Livelybrooks D, Jones P, et al. Causes of high crustal conductivity beneath the Iapetus suture zone in Great Britain. Geophys J Int , 1996, 124: 433-455. DOI:10.1111/gji.1996.124.issue-2 |
[5] | 陈小斌, 赵国泽, 马霄. 大地电磁三维模型的一维二维反演近似问题研究. 工程地球物理学报 , 2006, 3(1): 9–15. Chen X B, Zhao G Z, Ma X. Research on inversion of MT 3D model approximately by 1D, 2D inversion method. Chinese Journal of Engineering Geophysics (in Chinese) , 2006, 3(1): 9-15. |
[6] | Berdichevsky M N, Dmitriev V I. Magnetotellurics in the context of the theory of ill-posed problems. Tulsa:Society of Exploration Geophysicists, 2002 http://www.oalib.com/references/18987111 |
[7] | Berdichevsky M N, Dmitriev V I, Pozdnjakova E E. On two-dimensional interpretation of magnetotelluric soundings. Geophys J Int , 1998, 133: 585-606. DOI:10.1046/j.1365-246X.1998.01333.x |
[8] | Berdichevsky M N, Dmitriev V I, Kuznetsov V A. Bimodal two-dimensional interpretation of magnetotelluric sounding. Physics of the Solid Earth , 1996, 31(10): 821-837. |
[9] | Ledo J. 2-D versus 3-D magnetotelluric data interpretation. Surv Geophys , 2006, 27: 511-543. |
[10] | 肖骑彬.大地电磁测深反演方法对比与可视化及应用研究:中国科学院地质与地球物理研究所, 2004 Xiao Q B. Application study of magnetotelluric inversion method comparison and visualization. Beijing:Institute of Geology and Geophysics, Chinese Academy of Sciences, 2004 http://www.irgrid.ac.cn/handle/1471x/174762 |
[11] | Mackie R L, Madden T R, Wannamaker P E. Three-dimensional magnetotelluric modeling using difference equations--Theory and comparisons to integral equation solutions. Geophysics , 1993, 58(2): 215-226. DOI:10.1190/1.1443407 |
[12] | 陈小斌, 张翔, 胡文宝. 有限元直接迭代算法在MT二维正演计算中的应用. 石油地球物理勘探 , 2000, 35(4): 487–496. Chen X B, Zhang X, Hu W B. Application of finite-element direct iteration algorithm to MT 2-D forward computation. Oil Geophysical Prospecting (in Chinese) , 2000, 35(4): 487-496. |
[13] | 陈小斌, 赵国泽, 詹艳. MT资料处理与解释的Windows可视化集成系统. 石油地球物理勘探 , 2004, 39(Suppl.): 11–16. Chen X B, Zhao G Z, Zhan Y. A visual integrated windows system for MT data process and interpretation. Oil Geophysical Prospecting (in Chinese) , 2004, 39(Suppl.): 11-16. |
[14] | Rodi W, Mackie R L. Nonlinear conjugate gradients algorithm for 2-D magnetotelluric inversion. Geophysics , 2001, 66(1): 174-187. DOI:10.1190/1.1444893 |
[15] | 李墩柱, 黄清华, 陈小斌. 误差对大地电磁测深反演的影响. 地球物理学报 , 2009, 52(1): 268–274. Li D Z, Huang Q H, Chen X B. Error effects on magnetotelluric inversion. Chinese J. Geophys. (in Chinese) , 2009, 52(1): 268-274. |
[16] | deGroot-Hedlin C, Constable S. Occam's inversion to generate smooth, two-dimensional models from magnetotelluric data. Geophysics , 1990, 55(12): 1613-1624. DOI:10.1190/1.1442813 |