利用地磁场资料构建海域地磁场模型,在海底矿产资源勘探^[1]、地磁模型精确导航^[2]、水下目标探测^[3]、海洋环境污染^[4]等领域具有广泛的应用前景。

区域地磁场模型的构建方法主要包括球谐分析、球冠谐分析、等效源等方法^[5-6]。其中,球谐分析是当前全球和区域地磁场建模的主要手段,通过有限阶球谐函数的叠加来求解球谐系数,使地磁场数据达到最佳拟合^[7]。然而,由于球谐级数截断阶数的限制,该方法在空间分辨率方面存在局限性。例如,NGDC720模型尽管采用了720阶截断,但其空间分辨率仅为50 km,远不能满足精细探测的需求。此外,提高截断阶数虽然能提升分辨率,但会导致球谐系数的数量急剧增加,高低阶项系数相差悬殊,导致法方程矩阵高度病态,影响计算稳定性和解算精度^[8]。球冠模型在处理局部区域磁场时具有一定优势,但其建模效果依赖于球冠半角的选择,截断阶数的增大会导致级数收敛缓慢甚至失去物理意义,同时难以准确模拟磁场的径向变化,影响上下延拓的效果^[5]。

相比之下,等效源方法为区域地磁场建模提供了一种有效的替代方案。该方法通过在地下或地表设置等效源场源体,并求解其物性参数,使其能够拟合观测的磁场数据,进而重构磁场分布^[9]。依据高斯通量定理,理论上在满足位场不变的条件下,可以通过等效源重构地磁场模型,从而实现对区域地磁场的高精度建模。若等效源分布合理,重构磁场的精度可接近真实磁场分布。然而,该方法针对海域实际应用的特殊性仍然面临诸多挑战:首先,观测数据有限且离散,不合理的等效源布设可能会导致位场信号的“泄漏”,即尽管观测数据拟合精度较高,但未观测点的重构磁场误差较大^[10];其次,观测数据通常含有噪声,若等效源模型未能有效抑制噪声,可能会放大噪声影响,降低模型精度;另外,海域地质模型和观测条件的特殊性,也为等效源的构建增加了难度;最后,为了获得高精度和高分辨率的磁异常模型,需要构建大量等效源单元体,存储大规模系数矩阵,并求解大型线性方程组,计算效率和稳定性成为关键问题^[11]。因此,合理构造等效源布设方案,提高计算效率并增强模型稳定性,是等效源方法在区域地磁场建模中取得良好效果的关键。

本文围绕等效源方法在海域区域地磁场建模中的关键技术问题,探讨其理论基础、建模方法及优化策略,为构建高精度区域磁场模型及后续应用提供技术支撑。

等效源方法的实质是地球物理中的正演与反演问题^[12-13]。利用等效源方法构建区域地磁场模型可以分解为两步:第一步,首先通过拟合观测磁场数据建立等效源,即反演问题;第二步,利用等效源重构区域地磁场模型,即正演问题。具体可表示为:

式中:d为M×1维的观测数据向量;S为M×N阶的核函数矩阵;m为N×1维的模型参数向量。假设测点在地表规则分布,那么观测数据就可以表示为:

式中:ΔT_M表示第M个观测点的磁场值,M为观测点的个数。

假设等效源模型单元体在地下规则分布,那么等效源模型单元体的物性参数可以表示为:

式中:m_N表示第N个等效源模型单元体的物性参数,N为等效源模型单元体的个数。

求解式(1),可以表示为:

式中:m₁、m₂、…、m_N为未知参数,即等效源模型单元体的物性参数。基于获得的等效源参数,可以构建不同空间位置的三维地磁场模型。

理论上,若设置的等效场源与真实场源十分接近,将会得到高精度的重构与转换结果。由于反演计算具有多解性,尽管通过高精度拟合观测数据可以构建不同的等效场源,但不是所有的等效源都可以得到高精度的建模结果。一般使用均方误差RMS评价计算结果的精度:

式中:ΔT为理论数据;ΔT'为计算结果;N为数据量。但均方误差难以反映误差相对于理论值的大小,故采用理论值的振幅对均方误差进行归一化,得到归一化均方误差NRMS:

式中:ΔT_max为最大值理论;ΔT_min为最小值理论。

基于单个棱柱体模型分析等效源的等效性(模型一,图1a):模型尺寸为10 m×10 m×10 m,顶深为20 m,磁偏角0°,磁倾角50°,与地磁场方向一致;观测面位于0 m,网度为10 m×10 m。按照测网的网度设置了一层直立六面体等效源,单个等效源体的深度和尺寸与场源相同,等效源的磁化方向与场源相同。利用0 m观测面的数据反演等效源参数,进而计算上延5 m与下延5 m的结果。从图2a可以看出,重构得到的等效源模型能够完全重现场源特征,观测数据重构NRMS为9.39×10^-10;此时,上延NRMS为2.82×10^-10,下延NRMS为3.59×10^-9,两者与观测面重构结果基本处于同一误差水平。这表明当等效源设置合理时,其能够完全等效场源。

将等效源的深度设置上移1 m,此时等效源层位于实际场源体上方(图1b)。同样利用观测面的数据反演等效源参数,进而计算上延5 m与下延5 m的结果。此时,反演得到的等效源模型的中心位置磁化率降低,需要周围块体补充才能实现观测面数据的高精度重构(图2b4)。尽管观测数据重构NRMS仍可低至1.54×10^-9,但上延5 m与下延5 m误差明显增大(上延NRMS为2.80×10^-3,下延NRMS为1.97×10^-4)。这表明即使等效源位置设置存在微小偏差,等效源仍难以完全等效场源,进而导致延拓误差增大。因此,等效源如何设置尤为重要。

基于观测面以下介质的磁性特征,海平面以下大致分为3层:海水层、沉积层和地壳。海水层和沉积层的分界面为海底,沉积层和地壳的分界面为基底。其中,磁性成分主要分布在地壳中,海水层和沉积层基本无磁性。基于公开数据,建立如图3所示的海域磁性地质模型。其中,海底地形采用GEBCO-2020网格,空间分辨率为15″×15″;基底深度基于NGDC发布的GlobSedV3构建,空间分辨率为5'×5'。

从图3可以看出,距陆地区域大约100 km以内的大陆架及洋中岛屿附近的海水深度小于200 m,大陆坡海水深度介于200~2 000 m,洋底的海水深度介于2 000~6 000 m。而基底的深度主要由海水深度和沉积厚度控制,宏观上与海底相似,但在局部地区受火山活动和构造作用的影响,差异较大。根据不同的海底深度,设计了3个海域磁性简化模型(图4),分别为大陆架模型、大陆坡模型和大洋底模型。

海域地磁场三维建模的目标是利用船测或航磁观测数据构建海平面上、下一定范围内的三维地磁场模型,类似于位场数据的向上延拓与向下延拓。本文以海平面以上0~500 m高程地磁场建模为目标开展相关研究,并利用“上延”与“下延”的结果评价建模的精度。基于位场延拓理论分析,大陆坡与大洋底模型在250 m下延深度条件下,其延拓面与海底场源间距保持较大(延拓距离远小于源—观测面间距),使得向下延拓误差较小;然而大陆架模型因浅水地质特征导致同尺度向下延拓距离相对于源—观测面间距不可忽视,引发边界效应非线性放大问题,其最大相对误差显著增大,这成为浅海区高精度地磁场建模的关键技术瓶颈^[14]。

为了简化计算,一般采用简单的规则几何形体作为等效场源,例如:点元、面元和直立六面体。点元(偶极子)具有计算简单和设置方便的优点,但是偶极子为离散模型,对相邻场源中间区域难以约束。面元(偶层位)计算相对简单,但Li等^[15]指出偶层位与观测数据的实际物理特性无关。将均匀磁化直立六面体单元紧邻排列构成起伏等效层,其计算相对复杂,但更加符合物理基础^[10]。

在模型一的基础上,考虑一个深度更大的单个棱柱体模型(模型二),模型顶深为95 m,其他参数不变。分别利用偶极子等效源模型和直立六面体等效源模型进行延拓误差分析。从图5可以看出,直立六面体模型的延拓误差小于偶极子模型,尤其是上延误差。因此,采用直立六面体构建等效源模型更加符合实际地质情况,延拓结果具有更好的保幅性。

等效源的深度是等效源布设中的关键参数。图6为模型二等效源放置在不同深度时,延拓误差的变化曲线。从图中可以看出:等效源放置在不同深度时,对应的延拓精度有差异;当等效源设置深度与场源接近时,拟合误差与延拓误差最小,即存在最佳设置深度;且随着等效源深度与场源深度的偏差增大,误差有增大的趋势,即过浅或过深都会导致误差增大。

当等效源放置位置较浅时,尽管重构结果精度很高,但是延拓的误差较大。其原因是数学解不能保证浅部等效场源物性是光滑的。浅部的等效场源产生的异常为高频成分,随着高度的变化缩放很快,所以尽管对0 m平面测点数据拟合精度很高,但利用等效源计算加密网格(或未观测点)的数据仍存在高频误差(图7b)。且进行下延时,由于结果放大作用很快,导致延拓结果误差较大。当等效源放置位置较深时,由于无法拟合观测数据中的高频信息,无论是重构结果还是延拓结果,精度都较低。从物理角度来看,等效源位置过深时,无法高精度地拟合观测数据,进而无法实现高精度的延拓(图7a)。

因此,等效源深度参数的优化配置对位场延拓精度具有重要的影响,当其深度设置接近真实场源分布时,可得到高精度的重构和延拓结果,且当深度偏离场源实际位置时,延拓误差将呈非线性增长。特别地,等效源深度必须同时满足双重约束条件——既不能超过场源地质界面深度,又需大于5倍测点距,否则延拓结果中会产生高频噪声。

由于磁异常数据不仅与磁性场源的物性参数和几何参数有关系,还受到场源磁化方向的影响。利用模型二产生的磁异常数据,我们对比了不同磁化方向模型进行延拓的误差(图8)。其中,地磁场方向与等效源的磁化方向设置相同,均为50°。可以看出,当等效源磁化方向与真实场源磁化方向一致时(等效源顶深设置为95 m,图8a1),计算误差会达到最小;且随着等效源磁化方向与场源磁化方向的偏差增大,计算误差有增大的趋势。但需要指出的是,当等效源深度设置与场源深度一致时,合理设置等效源磁化方向的优势明显;但是当等效源深度设置与场源深度不一致时(等效源顶深设置为35 m,图8b2),磁化方向的影响变小。根据等效源深度设置结论,在实际工作中,为保证高频信息的拟合,等效源设置一般偏浅。因此,与等效源设置深度相比,磁化方向的影响明显更小。

考虑到自然界中存在的退磁与剩磁现象,会使地磁场方向与磁性场源真实的磁化方向不一致。将模型二中棱柱体场源模型的有效磁化方向的倾角设置为20°,地磁场仍为50°,开展等效源不同磁化方向延拓误差分析。结果表明,采用地磁场磁化方向和有效磁化方向设置等效源的延拓误差都不是最小的(图8b),最佳磁化方向介于有效磁化方向和地磁场磁化方向之间。

因此,与等效源深度设置相比,磁化方向的设置对延拓精度影响相对较小;当不存在剩磁影响时,等效源磁化方向设置为地磁场方向可以获得最小的延拓误差;当存在剩磁影响时,最佳磁化方向介于有效磁化方向和地磁场磁化方向之间。

当依据观测磁场反演获得等效源参数后,上延与下延具有不同的误差规律。通过分析影响上延与下延精度的因素,进而提高建模的精度。

建立顶深为295 m的单个棱柱体(模型三),尺寸为10 m×10 m×10 m,磁偏角0°,磁倾角50°。正演计算棱柱体模型在0 m平面的磁异常,观测网度为10 m×10 m。在设置等效源层时,每个测点下方放置一个直立六面体等效源。等效源的横向尺寸与测网相同,厚度为10 m,顶深215 m。通过反演等效源参数,进而计算上延0~200 m与下延0~200 m的磁异常。从图9a可以看出:①当延拓距离小于60 m时,上延与下延误差大致相当;②当延拓距离大于60 m后,上延误差趋于稳定,下延误差快速增大,即下延存在误差明显增大的转折点。因此,如何提高建模精度,尤其是降低下延误差是重要的问题。

针对模型三,我们测试了不同观测网度的结果。观测范围不变,观测网度分别为5 m×5 m和2 m×2 m。每个测点下方放置一个直立六面体等效源。

等效源的横向尺寸与测网相同,厚度均为10 m,顶深为215 m。由图9可以看出:不同网度的误差变化规律是相似的,且下延中误差增大转折点的位置并不受观测网度的影响。

基于模型三,通过改变场源的深度,我们测试了不同深度等效源的延拓结果,观测网度为5×5 m。测试中遵循本文提出的深度设置原则。从图10可以看出:不同等效源深度模型的计算结果具有相似特征,但随着等效源深度(顶深)的增大,转折点向后偏移,即能够保证下延精度的延拓距离与等效源的深度有关。通过对多个模型的结果进行线性拟合,在观测面附近一定范围内(0.23h,h为等效源的顶深)上延和下延精度大致相当。

针对模型三,我们测试了等效源不同厚度对延拓的影响,观测网度为5 m×5 m。从图11可以看出:与顶深相比,等效源的厚度对延拓精度的影响较小。但是,增大厚度能在一定程度上提高延拓精度。这表明更大的等效源厚度能够吸收更多长波信息,且在延拓过程中表现出更缓的衰减特性。

针对模型三,我们测试了不同扩边点数延拓的结果,观测网度为5 m×5 m,测点数20×20(图12)。可以看出:扩边可以显著提高延拓精度。精度提高体现在两个方面:第一是下延误差明显增大的拐点向后偏移,即上、下延精度大致相当的范围增大,扩边一倍后,范围大致可以从0.23h扩大至0.6h;第二是整体精度明显降低,按照5倍点距扩边后,相对误差小于1%。需要注意的是:从均方误差的角度看,下延误差大于上延误差,即上延的精度高于下延(图12a、b);从归一化均方误差的角度看,由于磁场上延存在明显衰减,下延误差可能小于上延误差(图12c、d)。

基于等效源方法的地磁场建模过程中,求取等效源的物性参数(式(4))是计算的主体,这是一项基于大量数据的最优化问题。其一,在等效源构建方面,核矩阵S内元素数据量为N²,存储数据会占用很大的内存空间,影响计算速度;其二,在应用角度方面,基于构建的等效源重建磁场模型时,调用全部的核函数和等效源参数也会导致计算速度的下降。传统方法将观测数据分割为若干个块,每个块单独进行等效源计算。以大陆架模型为例,将测区划分为25个区块,采用传统的分块方法进行试算。当延拓距离为50 m时,延拓结果与理论值相近;但当延拓200 m时,延拓结果已经完全失真(图13)。

为了解决简单分块方法带来的拼接问题,提高延拓精度,我们对分块计算方法进行了改进。如图14所示,采用窗口覆盖的方式进行分块计算:在进行反演计算时,为构建(延拓)区域A小范围的磁场数据,使用区域B中等范围的观测数据和区域C较大范围的等效源。计算完成之后,需要针对区域A保存区域C的等效源参数,并利用区域C的等效源参数实现区域A磁场的建模。区域C的等效源能够保证区域B观测数据的高精度重构。区域B观测数据能够保证区域A延拓数据的连续性。通过整体滑动,实现整个区域的建模。与整体算法相比,尽管窗口覆盖的方式会使很多数据重复参与计算,并增大了等效源参数的存储量,但是可以大大提高计算结果的精度。另外,该方法通过将整体数据剖分成多个小块,避免了整体核函数矩阵的计算与存储,在重构磁场时,也仅需要局部小窗口的等效源参与计算,提高了计算效率,尤其是在大规模数据的计算时,优势更为明显。

针对图4中的大陆架模型,采用N_c=7 Na,N_b=5 Na,进行窗口覆盖分块方法的试算。结果表明:下延50 m和上延50 m的归一化均方误差分别为0.16%和0.34%;下延200 m和上延200 m的归一化均方误差分别为0.42%和1.31%;这与整体算法的误差水平大致相当(下延0.33%;上延1.24%)。对比图13和图15可以看出,改进的分块方式的方式可以显著提高延拓的精度,减小窗口拼接导致的误差。

针对海域磁测特殊的地质和观测条件,结合等效源布设原则,采用具有一定厚度的直立六面体作为等效源。理论上,等效源的深度应小于真实场源的深度,且越接近真实场源,误差越小。考虑到海底地形资料的空间分辨率和精度均高于基底深度,因此沿着海底地形设置起伏的等效源层。以250 m飞行为例,当测线距抽稀到50 m间距时,在海底布设等效源层始终能够满足等效源深度大于5倍测点距的要求。

利用本文给出的等效源构建方案,针对图4中大陆架、大陆坡和大洋底3种模型开展试算。海底以下75 m存在两个尺寸均为100 m×100 m×10 m的棱柱体;深部存在1个尺寸为500 m×250 m×50 m的棱柱体,距海底顶深为1 155 m(图4);观测面位于海平面以上200 m。磁偏角0°,磁倾角50°。沿海底每个测点下方放置1个等效源体,单个等效源尺寸为20 m×20 m×1 000 m,等效源向4个方向分别扩边50个点。反演等效源参数,进而生成海平面以上0~500 m范围内的磁异常。另外,考虑到剩磁的影响,我们计算了等效源磁倾角为90°的误差,见表1所示。从计算结果可以看出,通过合理的场源布设,可以将误差控制在5%以内;其中在大陆坡区域,由于海底较浅,误差稍大;在大陆架和大洋底,随着海底深度的增加,建模精度较高。

1)通过对比分析发现,等效源几何参数与空间配置策略显著影响磁场模型精度。针对海域等效源地磁场建模的特殊性,详细分析了等效源模型的选择、模型参数的设置等问题,提出了地形跟随式直立六面体等效源的配置方案。

2)在算法实现层面,基于滑动窗口的覆盖式计算方案有效解决了大规模磁测数据的高精度处理瓶颈。实验结果表明,该方法既可保证拼接边界连续性,又能维持较好计算性能。