沈昭昂，彭 中，杨 波

(浙江大学 地球科学学院，浙江 杭州 310027)

1 引 言

大地电磁测深是深部地质结构探测和矿产资源勘查的重要方法，而海洋电磁法是近年来发展迅速的海底探测有效手段。在大地电磁观测中，由于相关噪声的影响，会对估算大地电磁张量阻抗产生干扰[1]。对海洋电磁探测而言，当观测区域邻近海洋时，洋流运动感生的电磁信号也会产生干扰。理论研究表明，海水运动产生的感应电磁场会影响后续的处理和解释结果[2]，而海水层也会对电磁场产生影响，主要表现为高频段的影响大于低频段,磁场的影响大于电场[3]。因此，当观测区域邻近海洋时，需要考虑海岸效应带来的影响[4]。

洋流运动主要包括潮汐、内波和湍流[5]。由地磁场中的海水运动引起的电磁场是海洋中大地电磁法或可控源电磁法的主要噪声源，特别是在水深1 km左右的沿海地区[5]。Kuvshinov等利用潮汐运动感生电磁信号来约束地球的电性结构[6]。Schnepf等通过分析数值模拟的结果，认为海底水平方向的磁场对岩石圈的导电层有着很高的敏感性，并将其归因于源区(海洋层)和海底介质之间的电流耦合[7]。Shimizu等尝试使用表面重力波产生的电磁场来探测海底下方的导电结构[8]，认为海洋长波产生的MT信号主要代表源的特性(海洋波或海啸)，不适用于探索地球深部地质结构。Grayver 等通过从卫星高度的M2潮汐磁数据的径向分量推断出一个新的全球电导率模型，证明了使用 M2潮汐场探索地球内部的可行性[9]。Zhang等利用全球3D电导率模型评估潮汐模型中的不确定性引起的OTEM(Ocean Tidal EM)数值预测误差，为岩石圈-软流圈系统内的良导结构提供了更多的约束[10]。

本文针对洋流运动产生的电磁信号是大地电磁正演中的噪声还是有效信号进行了研究。不考虑偶发海啸事件，正常海水运动感生磁场最大量级约为1 nT，然而在实际的大地电磁观测中，磁场的振幅最低只有0.001 nT[11]，这说明潮汐的运动可能会对海底和海岸带的大地电磁观测造成影响。本文选择了美国西海岸卡斯卡迪亚俯冲带为研究区，该海域的洋流为自北向南的加利福尼亚寒流，海水运动活跃，有着丰富的洋流数据用于海水感生电磁场的计算。该地区中新生代整体的构造环境以板块俯冲为主导，Juan de Fuca板块现今正沿北东向俯冲到北美板块下方，岩石圈电性结构高低电阻率对比强烈[12]，因而是研究海水运动感生电磁场源影响大地电磁响应的理想场所。此外，EarthScope在美国西海岸布设有大量MT(Magnetotelluric Method)观测台站，因此有着充足的大地电磁观测数据和大量的研究成果，这为建立正演的电性结构模型提供了充分的研究条件。本文建立了洋陆模型和卡斯卡迪亚俯冲带两种模型，并分别以电离层/磁层天然电磁场以及包含电离层/磁层天然电磁场和洋流运动的感生电磁场为场源，正演得到了两种模型的大地电磁响应。通过对比分析这两种场源情形下的大地电磁响应，获得了洋流运动感生电磁场对大地电磁响应的贡献。

2 洋流电磁效应模拟

2.1 海水运动产生的感应电磁场控制方程

考虑海水运动产生的感应电场问题，取时谐因子为eiωt，可以得到电场和磁场的控制方程：

(1)

(2)

(3)

对式(3)的等号两边同时求旋度，展开得到：

(4)

此处，外加场源的电流密度JS是洋流引起的，由地球磁场和洋流运动速度决定，表达为：

(5)

式(2)中，V为洋流的速度，m/s;σs为海水电导率，S/m;Bm为与时间无关的主要地磁场,T。式(4)和式(5)即为洋流电磁效应模拟的控制方程。有限体积法(Finite Volume Method,FVM)作为常用的求解偏微分方程的方法之一，通过将计算区域网格化剖分，在网格单元的网格点或者棱边处定义一个控制体积，将待求解的偏微分方程转化为对该控制体积求积分后的离散线性方程组。有限体积离散的核心是体积单元的设计，它的大小决定着计算的效率和精度。对每个待计算的电场分量构造一个体积单元，如图1所示。

图1 控制体积单元Ω及其表面边界Γ[13]Fig.1 Control volume elementΩ and its surface boundary Γ

绿色箭头标识的即待求解的电场分量，将式(5)向x,y,z轴进行投影，可以得到标量方程：

(6)

(7)

(8)

将海水速度以及地球主磁场的值与电场一样，定义在网格中心(图1)。式(6)～式(8)右端项离散化后在离散网格单元的内部每一个棱边上写出该电场分量所满足的方程，联立所有的方程，可以得到一个稀疏的复对称线性方程组：

(9)

2017年6月底，大同市市级和7县4区的实施方案已全部出台，方案的制定和出台为河长制的全面推行提供了行动指南。建立市、县、乡、村四级河长体系，将河长制拓展到了村一级。各级河长上岗履职，积极开展巡河治河工作。截至2017年底，大同市共竖立河长公示牌453块，将流域面积大于50 km2的106条河流全部纳入河长制体系，构建了层级清晰、责任明确、推进有力的河长制组织领导体系，形成了一级抓一级、层层抓落实的组织管理体系。

2.2 场源的处理

有限体积法求解海水运动的控制方程，得到的是海水以真实的运动速度切割地磁场产生的感生电磁场。大地电磁正演模拟得到的是阻抗、相位、电阻率等参数，为了计算的简化,通常假定第一层电场的幅值为1，最后得到的电磁场值为相对值。若要研究海水运动感生电磁场对大地电磁场观测数据的影响，可以将大地电磁正演得到的相对电磁场值与真实的电磁场幅值量级进行对比，将其转化为真实的量级，进而与海水感生的电磁场叠加形成新的大地电磁场观测值，最后计算其对应的阻抗、视电阻率及相位。对于电场相对值与真实值的转换，首先需要验证一个假设，即在正演中，当模型第一层输入的电场值为1和n时，最后分别得到电场E1和E2，那么有E2=nE1。以x极化下大地电磁三维正演电场的求解为例，在大地电磁正演中当第一层输入的电场值为1(上边界条件u=1),最后一层边界条件为u=0，方程的右端项即为首项为1的矩阵b1；当第一层输入的电场值为n(上边界条件u=n),最后一层边界条件为u=0，方程的右端项即为首项为n的矩阵b2。在两种边界条件情形下，电场的代数方程式分别为AE1=b1,AE2=b2，因为b2=nb1，可以证明电场E1和E2满足关系E2=nE1。

本文主要根据全球电磁场强度平均振幅谱特征图(图2)来进行电磁场真实值与相对值的转换。从图2中可以看出，大地电磁场的信号很微弱，电场振幅最低只有0.01 mV/km，磁场的振幅最低仅为0.001 nT。对于如此微弱的信号，很容易被环境噪声掩盖,导致信号失真。根据图2可以将大地电磁正演模拟得到的电磁场值转化为相应周期的真实幅值的量级，从而可以与真实的海水运动产生的电磁响应叠加，来分析海水的运动产生的电磁效应是否会对大地电磁微弱的观测信号产生影响。

图2 全球电磁场强度平均振幅谱特征[11]Fig.2 The characteristics of the average amplitude spectrum of the global electromagnetic field intensity

电磁场含不同周期的信号，当计算某个观测周期下海水运动的感生电磁场值时，需要找到对应周期的大地电磁观测的真实幅值。当选择周期为1 000 s的信号时，电场的平均振幅为2 mV/km；当选择周期为1 s的信号时，电场的平均振幅为0.2 mV/km。

2.3 正确性检验

验证数值解的正确性需要与解析解进行对比。对于电偶极子源，Cox给出了一维电导率模型的解析解[14]。通过计算半空间海水模型(电导率3.3 S/m)和一维层状模型(3.3 S/m，0.1 S/m，0.2 S/m)的有限体积解，与解析解进行对比验证。解析解的计算采用加州大学圣地亚哥分校Scripps海洋研究所(Scripps Institution of Oceanography,SIO)海洋电磁实验室所开发的Dipole 1D进行计算[15]。对于电偶极子有限体积法数值模拟，需要考虑源的奇异性问题，本文使用了伪delta函数[16]的方法计算可控源电磁响应。通过一个小的区域代替一个点源，消除了源的奇异性。伪delta函数的表达如式(10)所示。

(10)

式(10)中，x0表示场源点，τ为常数，决定了伪delta函数的幅值大小和加载节点的数量。三维情况下偶极子源可以表示为：

Js=I0·dl·δs(r-r0)

(11)

其中，

δs(r-r0)=δs(x-x0)δs(y-y0)δs(z-z0)

(12)

式(11)中，I0为电偶极子源电流大小,A;dl为电偶极子长度,m。一维层状模型中，电偶极子设置在(0, 0, 100)，偶极子长度为100 m，电偶极距为1 A。m为发射周期为1 s，发射方向为x方向。接收站分布于x=0 m,y向从-200～1 000 m,z=0 m，接收点距为100 m。设计的网格为Nx=33,Ny=33,Nz=48，网格节点数是52 272，自由度数为156 816。计算区域水平方向大小为2 000 m×2 000 m，网格间距100 m。然后逐渐向两边扩展。垂向网格在地表加密，随着深度增加网格逐渐变大。为了提高对于电偶极子模拟的精确度，在电偶极子所在的位置进行网格加密。电导率结构为包含海水(3.3 S/m)、沉积地层(0.1 S/m)及岩石圈(0.02 S/m)，水平电场分量Ex的幅值如图3所示。

图3 一维层状模型的电偶极子水平电场幅值曲线3D FVM解与解析解对比Fig.3 Comparison of the 3D FVM solution and the analytical solution of the electric dipole horizontal electric field amplitude curve of the one-dimensional layered model注：图中，y轴左侧坐标为电场对数域的幅值，右侧坐标为有限体积解与解析的相对误差。蓝色实线为一维解析解，蓝色虚线为有限体积解，洋红色虚线为两者的百分相对误差

从图3中可以看到，由于加载了伪delta函数来模拟电偶极子，因此在源区电场响应的拟合结果较好，在远场有限体积解与解析解一致性很好，相对误差小于2 %，可以满足后续对比研究的需要。

3 洋流电磁效应的数值模拟

3.1 研究区数据条件

在研究区海水感生电场的模拟中，磁场数据采用国际地磁参考场(International Geomagnetic Reference Field, IGRF)模型[17]，海水速度模型从HYCOM(HYbrid Coordinate Ocean Model)获取[18]。其中，IGRF模型计算了从1900年至今的地球主磁场。根据IGRF-13模型，研究区域俄勒冈州海岸带的磁场及其变化如图4所示，由于海水仅在0～5 km有可测的速度，因此研究区磁场分量的垂向分布只在5 km以浅，且磁场在5 km深度范围内变化很小，2 m深度和5 km深度的磁场值非常接近，大约相差几十个nT，这对于主磁场的量级来说，影响很小。

图4 美国西海岸磁场Fig.4 US west coast magnetic field

HYCOM模型提供了全球HYCOM+NCODA的近实时海洋预测系统。本文选取GOFS 3.1: 41-layer HYCOM+ NCODA Global 1/12° Analysis模型数据，按照GLBy0.08的网格分布，即在经度上精度为0.04°，在纬度上精度为0.08°，海水的周期是3小时，数据格式采用NetCDF。研究区域不同深度的海水速度切片如图5所示，海水速度结构是二维的，即在水平方向上分布，图中箭头表示速度矢量，由东向和北向的速度合成。在海水深度100 m以浅，有着丰富的洋流活动，从100 m开始海水的速度大幅降低，到400 m只在0～0.1 m/s左右，可以说其对物理场的模拟影响非常小。从图5(f)中可进一步看到，5 km深度处的海水几乎不存在区域尺度上的速度结构。

图5 研究区不同深度的海水速度Fig.5 Seawater velocity at different depths in the study area

3.2 简单洋陆模型

在美国西海岸，岩石圈的电阻率高达几千甚至上万欧姆米，而海水的电阻率仅为0.3 Ω·m。对于实际的离散线性方程组来说，左端项的病态程度会很大，因此需要对该电导率模型进行简化。参考全球和区域电导率分布[19]，建立简化洋陆模型(图6)。海洋一侧组成岩石圈的三层电阻率为0.3 Ω·m—50 Ω·m—100 Ω·m，陆地一侧组成岩石圈的三层电阻率为100Ω·m—50 Ω·m—100 Ω·m，下伏软流圈电阻率为10 Ω·m，地幔转换带电阻率为25 Ω·m，地幔转换带底界面660 km以下电阻率为3 Ω·m。

图6 简单洋陆模型及网格划分Fig.6 Simple ocean-continent model and meshing

利用有限体积法求解大地电磁三维正演问题，网格剖分的大小疏密都会对计算的结果产生影响。无穷远边界处的近似条件会加大对计算结果的影响，而边界距离测点越近，影响就会越大。在有限计算能力和尽可能保证计算效率的情况下，通常将计算区域划分为两个区域：目标区域和网格外延区域。目标区域一般是包含测点的地质体存在的核心区域，而网格外延区域则是越靠近边界越稀疏。这样的网格划分是符合电磁场的物理衰减规律的。外延区域的网格优化可以避免很多不必要的计算，在保证精度的同时节省计算量。网格剖分如图6(b)和图6(c)所示。网格的大小和稀疏不仅会影响求解方程的收敛性以及计算速度，还会影响正演的精度。当电磁场以平面波形式传播时，正演的误差与水平方向的网格有关[20]。因此对模型进行有限体积网格的剖分时，核心区域内(200 km×200 km)网格均匀，纵向则在浅部加密，深部稀疏，这样可以得到较为可靠的地表的观测值。对于外延区域，每一次扩边设计七个单元网格为一组，按照给定的系数自动进行扩大，越靠近边界，网格越大。

根据磁场模型和海水速度模型，通过有限体积法可以计算研究区域洋流运动产生的感生电磁场。计算的周期范围为100～1 000 s。选取100 s和 1 000 s周期，截取400 m深度的水平切片，它们在x极化和y极化下的MT水平电场如图7和图8所示。

图7 两种情况下，观测周期为1 000 s的水平电场Fig.7 The horizontal electric field with the observation period of 1 000 s in two cases

图8 两种情况下，观测周期为100 s的水平电场Fig.8 The horizontal electric field with the observation period of 100 s in two cases

图7和图8中可以看到，两个不同深度的剖面都表明在长周期信号下，水平电场在东西方向呈条状分布，变化缓慢。在观测周期为100 s时，水平电场仅在海岸线附近呈现较大的变化，而在海水或者陆地区域变化较小(图8a,图8b)。图8(b)和图8(d)的等值线表明，海水感生电磁场对y极化的水平电场比x极化的水平电场影响更大，对y极化电场的影响集中在-50 km

设计S1测线沿东西走向(-100～100 km)，每隔10 km布设一个大地电磁观测台站，其中海洋区域的台站设计在海底(2.5 km)处。不同周期的大地电磁响应结果如图9所示。

图9 两种情况下，不同周期的视电阻率及相对误差变化曲线Fig.9 Apparent resistivity and relative error curves of different periods in two cases注：下标ModEM表示不考虑海水运动情况下的计算结果，相应的分量用散点符号表示；下标Ocean表示考虑了海水运动情况下的计算结果，相应的分量用实线表示; 其中蓝色表示xy分量，红色表示yx分量。下标error表示两种情况下对应分量的百分相对误差并用虚线表示，浅蓝色代表xy分量的相对误差，紫色表示yx分量的相对误差。

图9和图10分别为东西向剖面的大地电磁观测的视电阻率以及相位变化，右侧坐标是考虑和不考虑海水运动两种情况下大地电磁响应之间的百分相对误差。在陆地一侧，相对误差较小(<0.1 %)，海水运动的电磁效应对大地电磁观测几乎不产生影响。在海洋一侧-50～-40 km处，洋流对相位影响较大，这或是由于该处的相位本身处于一个较小值以及该处的洋流速度较大，对电磁场的相位产生了较大的影响，但整体相位的趋势和未加入海洋源的相位保持一致。洋流电磁效应对视电阻率的影响范围比对相位的影响范围要更大，横跨海岸线到-50 km以远，在-40 km和-50 km处相对误差高达约50 %。在低于1 000 s时，海水运动的电磁效应会使得视电阻率增大；而当周期为1 000 s时，海水运动的电磁效应使得视电阻率减小。当观测信号在200 s附近时，海水运动产生的电磁场对于大地电磁响应的影响最大，这可能与该频段天然电磁场信号强度较低有关。此外，相位相比于视电阻率对海水运动产生的电磁场更为灵敏，在海水最为活跃的区域(-50～-40 km)，在较短的周期就出现相位异常变化，但出现视电阻率异常变化的周期要长。综上所述，电磁效应对陆边区域的大地电磁观测影响较小，但在海水区域50 km范围内的大地电磁观测的影响较大。

图10 两种情况下，不同周期的相位及相对误差变化曲线Fig.10 Phase and relative error change curve of different periods in two cases注：下标ModEM表示不考虑海水运动情况下的计算结果，相应的分量用散点符号表示；下标Ocean表示考虑了海水运动情况下的计算结果，相应的分量用实线表示; 其中蓝色表示xy分量，红色表示yx分量。下标error表示两种情况下对应分量的百分相对误差并用虚线表示，浅蓝色代表xy分量的相对误差，紫色表示yx分量的相对误差。

3.3 卡斯卡迪亚俯冲带模型

Egbert等通过对大地电磁数据进行反演得到了美国Cascadia俯冲带的三维电阻率模型[21]。该模型成功恢复出俯冲的大洋岩石圈，能够较好地反映美国西海岸岩石圈电性结构。本文截取该电阻率模型的核心区域作为正演模型，设计了一条横跨海洋和陆地的剖面，在考虑海水运动感应电场的情况下进一步进行三维正演模拟，以探究海水运动感生电磁场对于美国西海岸大地电磁观测的影响(图11)。其中，图11(a)为美国西海岸地形，图中红色倒三角为设置的测线S1及台站分布，海上的台站设置在海底，随海底地形变化，陆地的台站设置在地表，右下角红色矩形指示的是研究区域。图11(b)～图11(e)为具体的正演模型及网格划分，网格的核心区域在水平面设置200 km×200 km，深度到800 km。图11(c)中沿S1测线每隔10 km布设一个台站，其中Js表示海水运动产生的感生电磁源。

图11 美国西海岸地形及正演模型设计Fig.11 U.S. west coast terrain and forward modeling design

选取模型采用迭代125次的反演结果，因此模型具有很好的光滑性，在低频信号下求解控制方程具有很好的收敛性，计算的周期范围为20 ～10 000 s，选取20 s和10 000 s周期，截取400 m深度的水平切片，它们在x极化和y极化下的水平电场如图12和图13所示。

图12 两种情况下，观测周期分别为10 000 s的水平电场Fig.12 The horizontal electric field with the observation period of 10 000 s in two cases

图13 两种情况下，观测周期分别为20 s的水平电场Fig.13 The horizontal electric field with the observation period of 20 s in two cases

从图12(b)、图12(d)、图13(b)、图13(d)中可以看到，水平电场Ey在400 m深度的等值线切片上，很好地分辨了海洋和陆地的分界线。短周期信号(20 s)在x极化下水平电场在东西向的衰减要比长周期信号(10 000 s)快，短周期信号下水平电场在-25 km处大幅度减小，这可能与俯冲带模型在该处俯冲板块上方存在南北走向的低阻体(图11b)有关。当观测周期大于1 000 s时，考虑和不考虑海水运动两种情况下的电场没有太大的差异，海水运动产生的扰动对于大地电磁观测的影响很小(图12)；但当观测周期小于1 000 s时，x极化的电场在50 s观测周期下，在海洋一侧-10～-100 km区域，存在多处电场的高值(图13c, 图13d)，这是由于海水运动感生电磁源带来的电场异常，在这些电场异常的区域，海水运动活跃、运动速度大，这些异常可能会对海底观测台站产生较明显的干扰。

设计S1测线(图11a)沿东西走向(-100～100 km)，每隔10 km设计一个大地电磁观测台站，分析海水产生的电磁场对大地电磁观测的影响，其中海洋区域的台站设计在海底，根据地形变化，最深可达1.5 km。不同周期的大地电磁响应结果如图14所示。

图14 两种情况下，不同周期的视电阻率及相对误差变化曲线Fig.14 Apparent resistivity and relative error curves of different periods in two cases注：下标ModEM表示不考虑海水运动情况下的计算结果，相应的分量用散点符号表示；下标Ocean表示考虑了海水运动情况下的计算结果，相应的分量用实线表示; 其中蓝色表示xy分量，红色表示yx分量。下标error表示两种情况下对应分量的百分相对误差并用虚线表示，浅蓝色代表xy分量的相对误差，紫色表示yx分量的相对误差。

图14和图15分别为50 ～10 000 s观测周期下，东西向剖面的大地电磁观测的视电阻率以及相位变化，右侧坐标是考虑和不考虑海水运动两种情况下大地电磁响应之间的百分相对误差。当观测周期为1 000～10 000 s的低频信号时，海水运动对陆地观测的影响很小，相对误差接近于零，而在海水区域尤其是海水运动速度大的地区，海底大地电磁观测的相对误差最大会达到20 %。当观测周期小于1 000 s时，陆地一侧的大地电磁观测受海水运动影响很小，相位不发生偏转，但是在海水一侧(-70～-10 km)电场的相位发生了转大的偏转，视电阻率出现了大幅度的降低，这比没有海水运动感生电磁源的大地电磁模拟得到的观测结果要低一个数量级，这可能与高频信号下电场的幅值低有关。因此，海底大地电磁观测信号可能会受到洋流的干扰，使得观测信号的可靠性大幅下降。

图15 两种情况下，不同周期的相位及相对误差变化曲线Fig.15 Phase and relative error curves of different periods in two cases注：下标ModEM表示不考虑海水运动情况下的计算结果，相应的分量用散点符号表示；下标Ocean表示考虑了海水运动情况下的计算结果，相应的分量用实线表示; 其中蓝色表示xy分量，红色表示yx分量。下标error表示两种情况下对应分量的百分相对误差并用虚线表示，浅蓝色代表xy分量的相对误差，紫色表示yx分量的相对误差。

4 结 论

本文通过海水运动的大地电磁正演模拟研究了洋陆边界简化地电模型和卡斯卡迪亚俯冲带模型的海水感生电磁场对大地电磁正演的影响。在洋陆边界简化地电模型中，海水运动产生的电磁场对于陆基台站观测几乎不产生影响，相对误差在0.1 %以下。而在海洋海水运动速度大的区域，其感生的电磁场对于大地电磁的观测不可忽略，相对误差可高达50 %。对于观测信号而言，200 s附近的信号观测到的大地电磁场受到的影响最大，这可能既与海水运动的速度有关，也与该频段天然电磁场源信号强度低有关。在模拟信号的周期范围内，当观测信号从低频变为高频时，受到的海水运动感生电磁场的影响逐渐增大。在卡斯卡迪亚俯冲带模型中，海水运动产生的电磁场对于陆基台站观测几乎不产生影响，相对误差接近于零。在观测频段为1 000～10 000 s时，在海水运动活跃的区域，海底大地电磁观测的相对误差大都在10 %以下。当观测周期小于1 000 s时，在靠近海岸的海水运动活跃区域，电场的相位发生较大偏转，视电阻率的值降低了约一个数量级，海底台站观测信号受到海水感生电磁场的影响强烈。因此，本文认为洋流的电磁效应对于海底观测台站会产生较大的干扰，而对陆基台站的观测所带来的影响可以不予考虑。

致谢

感谢浙江大学徐义贤教授对本文提出的的宝贵建议;感谢美国加州大学圣地亚哥分校海洋学院海洋电磁实验室提供的Dipole 1D程序以及弗吉尼亚理工学院的 Chester Weiss在程序编写等方面提供的帮助。