CN118033764B - 一种多尺度多方法地质体电阻率三维成像方法及系统 - Google Patents

Abstract

本申请提供了一种多尺度多方法地质体电阻率三维成像方法及系统，属于地球物理电磁法勘探领域，方法包括：对地下空间使用任意剖分方法进行离散后获取反演网格，构建初始反演电阻率模型；对不同类型的观测数据进行分组；将当前反演电阻率模型转换至为各个分组设计的局部正演网格，采用相应的数值求解方法分组分别进行正演，将理论正演响应数据与观测数据进行求差，获取数据拟合差；若数据拟合差大于预设拟合差，则计算缩减雅各比矩阵；将缩减雅各比矩阵映射到反演网格，再使用梯度类方法迭代更新反演模型。本申请可以实现同步并行处理不同类型数据，使用不同数值求解方法处理不同尺度的观测数据。

Description

一种多尺度多方法地质体电阻率三维成像方法及系统

技术领域

本申请属于地球物理电磁法勘探领域，更具体地，涉及一种多尺度多方法地质体电阻率三维成像方法及系统。

背景技术

现有的地质体三维电阻率反演方法存在以下问题：

一、难以适应多方法联合成像的需求：电阻率三维成像方法中获取电磁数据时采用的观测手段包括地球物理电法和电磁法；电磁法又包含各类细分方法，且都对地下空间电阻率分布敏感，比如电磁法分为直流电磁方法、频率域电磁方法和时间域电磁方法；而每一种电磁法又包括不同的激发源类型、不同的接收观测场类型、不同的观测模式和架设条件。例如：为时间域电磁法设计的反演程序需要在时间域求解扩散场形式的麦克斯韦方程组，而为直流电磁方法设计的反演程序需要求解稳恒状态下的泊松方程。虽然单独观测手段的三维反演已经日趋成熟，但是两种观测手段联合反演的方法难以处理任意多方法联合电阻率成像的需求。主要问题在于现有方法仅能处理一种或两种事先固定的数据类型，一旦需要处理新的数据类型，就必须重新编制程序且能接受的数据类型必须固定不变。

二、难以适应多尺度反演成像的需求：不同的电磁法具有不同的探测深度和成像尺度；例如，直流电磁方法的成像尺度和探测深度为几米至几十米；瞬变电磁方法的成像尺度为几十米至几百米，天然源大地电磁方法的成像尺度为几千米至几十千米；即使在同一种方法内部，也同时存在多个不同的成像尺度。传统反演方法反演多个成像尺度时需要使用繁杂的网格剖分技术同时兼顾小尺度网格和大尺度网格，这就造成了巨量的网格单元数，而大量的网格显著降低求解速度，甚至导致现有计算资源不可解。

三、难以满足不同数值混合求解的需求：例如，在地形起伏或地质目标体结构复杂的场合下，基于任意四面体网格的有限元法精度更高；而对于地形平坦且背景简单目标体的场合，积分方程法甚至半解析解计算效率更高。一个地质勘探项目可能跨越不同的地质背景条件，所以一个大的区域模型最好针对不同位置和地质条件使用不同的数值求解技术，然而，现有方法提前预置了网格离散化方法和数值求解方法，不能很好满足实际问题中地质条件复杂多变，因地制宜的需求。

发明内容

针对现有技术的缺陷，本申请的目的在于提供一种多尺度多方法地质体电阻率三维成像方法及系统，旨在解决现有地质体三维电阻率反演方法无法使用不同数值混合求解，导致难以适应多观测方法多尺度联合成像的问题。

为实现上述目的，第一方面，本申请提供了一种多尺度多方法地质体电阻率三维成像方法，包括以下步骤：

S1：对地下空间使用任意剖分方法进行离散后获取反演网格，并在反演网格中填充电阻率值，构建初始反演电阻率模型作为当前反演电阻率模型；

S2：对不同类型的观测数据进行分组，以便分组后的每组观测数据对应一个正演网格和正演模型，使用相同数值求解方法进行正演；

S3：将当前反演电阻率模型映射至正演网格，再根据每组观测数据对应的数值求解方法进行正演，计算每组观测数据对应的理论正演响应数据，并将理论正演响应数据与观测数据进行求差，获取数据拟合差；

S4：若数据拟合差小于等于预设拟合差，则判定当前反演电阻率模型能够用于地质体电阻率三维成像；否则转至S5；

S5：基于缩减雅各比矩阵和正演网格到反演网格之间的映射关系，完成定义在反演网格上的完整雅各比矩阵与任意向量的乘积操作，求取当前反演电阻率模型的总更新向量；其中，所述缩减雅各比矩阵的行数是正演模型对应分组内观测数据的个数，列数为分组对应正演网格的单元个数；其中，所述任意向量中包含有数据拟合差信息；

S6：将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型；转至S3，直至数据拟合差小于等于预设拟合差或者迭代次数达到最大迭代次数，判定当前反演电阻率模型能够用于地质体电阻率三维成像。

进一步优选地，步骤S2具体包括以下步骤：

根据获取观测数据采取的不同观测手段，将所述观测数据分为不同类型的观测数据集；

再将每一组观测数据集内的观测数据根据位置的相似性、地质结构相似性和电磁感应尺度的相似性进一步划分；

根据划分后的每组观测数据确定对应的正演网格和正演模型以及不同的数值求解方法。

进一步优选地，将当前反演电阻率模型映射至正演网格的方法为：

将每个正演网格单元与反演网格单元相交的体积作为权重，构建映射矩阵；

将反演网格单元内电阻率按照映射矩阵对正演网格单元内的平均电阻率进行贡献。

进一步优选地，S5中求取当前反演电阻率模型的总更新量向量的方法为：

S5.1：求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

S5.2：采用映射矩阵与各分组缩减雅各比矩阵相乘，获取各分组对应的完整雅各比矩阵；

S5.3：多次迭代搜索并综合各个分组完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更新向量；

或

S5.1：求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

S5.2：多次迭代搜索任意向量，并将任意向量与映射矩阵相乘，得到长度与正演模型单元个数相同的缩减向量；

S5.3：采用缩减向量与各正演模型对应的缩减雅各比矩阵相乘，获得各分组对应完整雅各比矩阵与任意向量相乘的结果；

S5.4：综合各分组对应完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更新向量。

进一步优选地，当采集的观测数据来自平原地区且地层成层性超过预设成层性时，正演模型为层状模型，数值求解方法为解析解或半解析解求解方法；

当采集的观测数据来自山地且山体有异常体时，反演网格为非规则网格剖分，数值求解方法为有限元求解方法；

当采集的观测数据来自平原地区且地下存在局部不均匀体时，反演网格为规则网格剖分，数值求解方法为有限差分、有限体积或积分方程方法。

进一步优选地，当反演网格和所述正演网格不使用同一坐标系时，获取正演模型与当前反演电阻率模型的电阻率对应关系时进行坐标转换。

第二方面，本申请提供了一种多尺度多方法地质体电阻率三维成像系统，包括：

初始反演电阻率模型构建模块，用于对地下空间使用任意剖分方法进行离散后获取反演网格，并在反演网格中填充电阻率值，构建初始反演电阻率模型作为当前反演电阻率模型；

观测数据分组模块，用于对不同类型的观测数据进行分组，以便分组后的每组观测数据对应一个正演网格和正演模型，使用相同数值求解方法进行正演；

数据拟合差获取模块，用于将当前反演电阻率模型映射至正演网格，再根据每组观测数据对应的观测手段，采用相应的数值求解方法进行正演，计算每组观测数据对应的理论正演响应数据，并将理论正演响应数据与观测数据进行求差，获取数据拟合差；

反演电阻率模型判定模块，用于判断数据拟合差是否小于等于预设拟合差；若数据拟合差小于预设拟合差，则判定当前反演电阻率模型能够用于地质体电阻率三维成像；同时用于判断当前迭代次数是否达到最大迭代次数，若达到最大迭代次数则判定当前反演电阻率模型能够用于地质体电阻率三维成像；

反演模型的更新模块，用于若数据拟合差大于预设拟合差，则计算各正演模型对应的缩减雅各比矩阵；基于缩减雅各比矩阵和正演网格到反演网格之间的映射关系，完成定义在反演网格上的完整雅各比矩阵与任意向量的乘积操作，用以求取当前反演电阻率模型的总更新量向量；将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型；

其中，任意向量中包含有数据拟合差信息；缩减雅各比矩阵的行数是正演模型对应分组内观测数据的个数，列数为分组对应正演网格的单元个数。

进一步优选地，观测数据分组模块包括：基于观测手段的分组单元、基于相似性的分组单元和正演方法确定单元；

基于观测手段的分组单元用于根据获取观测数据采取的观测手段不同，将所述观测数据分为不同类型的观测数据集；

基于相似性的分组单元用于将每一组观测数据集内的观测数据根据位置的相似性、地质结构相似性和电磁感应尺度的相似性进一步划分；

正演方法确定单元用于根据划分后的每组观测数据确定对应的一个正演网格、正演模型以及一种数值求解方法。

进一步优选地，数据拟合差获取模块中将当前反演电阻率模型映射至正演网格的方法为：

将每个正演网格单元与反演网格单元相交的体积作为权重，构建映射矩阵；

将反演网格单元内电阻率按照映射矩阵对正演网格单元内的平均电阻率进行贡献。

进一步优选地，反演模型的更新模块包括：缩减雅各比矩阵求解单元、分组完整雅各比矩阵求解单元、更新向量求解单元和相加单元；

缩减雅各比矩阵求解单元用于求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

分组完整雅各比矩阵求解单元用于采用映射矩阵与各分组缩减雅各比矩阵相乘，获取各分组对应的完整雅各比矩阵；

更新向量求解单元用于多次迭代搜索并综合各个分组完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更新向量；

相加单元用于将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型；

或

所述反演模型的更新模块包括：缩减雅各比矩阵求解单元、缩减向量求算单元、更新向量求解单元和相加单元；

所述缩减雅各比矩阵求解单元用于求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

所述缩减向量求算单元用于多次迭代搜索任意向量，并将任意向量与映射矩阵相乘，得到长度与正演模型单元个数相同的缩减向量；

所述更新向量求解单元用于采用缩减向量与各正演模型对应的缩减雅各比矩阵相乘，获得并综合各分组对应完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更细向量；

相加单元用于将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型。

第三方面，本申请提供一种电子设备，包括：至少一个存储器，用于存储程序；至少一个处理器，用于执行存储器存储的程序，当存储器存储的程序被执行时，处理器用于执行第一方面或第一方面的任一种可能的实现方式所描述的方法。

第四方面，本申请提供一种计算机可读存储介质，计算机可读存储介质存储有计算机程序，当计算机程序在处理器上运行时，使得处理器执行第一方面或第一方面的任一种可能的实现方式所描述的方法。

第五方面，本申请提供一种计算机程序产品，当计算机程序产品在处理器上运行时，使得处理器执行第一方面或第一方面的任一种可能的实现方式所描述的方法。

可以理解的是，上述第二方面至第五方面的有益效果可以参见上述第一方面中的相关描述，在此不再赘述。

总体而言，通过本申请所构思的以上技术方案与现有技术相比，具有以下有益效果：

本申请提供了一种多尺度多方法地质体电阻率三维成像方法，其中，对来自多个不同观测数据集进行分组，便于同组的观测数据使用相同的方式求解，而组和组之间可以不同，因此可以实现并行处理不同类型数据，使用不同数值求解方法处理不同尺度的观测数据。

本申请提出了一种多尺度多方法地质体电阻率三维成像方法，分别为正演模型和反演模型设置反演网格和正演网格；反演网格和反演模型仅用于统筹不同观测数据以及反演结果，不用于数值计算，因此反演网格可以任意复杂（细网格、大范围和任意形状单元），不增加求解困难；具体数值求解由多个正演网格和正演模型完成，因此各个数据分组可以根据其具体条件和需求采用不同的最优方法，整体上可以达到性能的最优，不会陷入一个方法难以适应多种需求的困境。如果数据分组足够细，则正演网格单元数量非常小，因此有利于使用直接求解法计算麦克斯韦方程离散得到的线性方程组；对于小型线性方程组而言，直接求解比迭代解法计算速度快，计算精度高，内存消耗适中。

本申请提出了一种多尺度多方法地质体电阻率三维成像方法，使用一个线性映射矩阵来处理从反演模型到正演模型的映射，以及从缩减雅各比矩阵到完整雅各比矩阵的映射；这两个映射互为转置关系，仅需要保存一个映射矩阵即可，得益于这一映射矩阵，反演过程中计算量最大最复杂的两个步骤（正演计算和雅各比矩阵与反演模型中的电阻率向量相乘）都可以在各个分组内完成，分组间互相独立便于并行求解，具有良好的可扩展性。

现有的电阻率模型的获取正反演共用的是一套网格的思路，本申请按照功能区分为反演网格和正演网格，且正演网格和正演模型可以根据数据分组为多个，互相独立；反演网格包含完整区域，有多个不同尺度、兼容各种地形地质和观测模式，每个正演网格仅负责一小部分数据、单一的尺度、单一的观测模式和某种特定的数值方法。

附图说明

图1是本申请实施例提供的反演方法流程图；

图2是本申请实施例提供的数据分组与反演、正演网格示意图；

图3是本申请实施例提供的多尺度多方法地质体电阻率三维成像系统。

具体实施方式

为了使本申请的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本申请进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本申请，并不用于限定本申请。

在本申请实施例中，“示例性的”或者“例如”等词用于表示作例子、例证或说明。本申请实施例中被描述为“示例性的”或者“例如”的任何实施例或设计方案不应被解释为比其它实施例或设计方案更优选或更具优势。确切而言，使用“示例性的”或者“例如”等词旨在以具体方式呈现相关概念。

首先，对本申请实施例中涉及的技术术语进行介绍。

观测数据：观测的地球物理场数据，例如：电场、磁场、或其简单导出量（如磁场的时间导数、时间序列傅里叶变换后的系数、电场与磁场的比值（阻抗）等）；狭义的数据指将所有待解释的地球物理场数值按照某种顺序排列成一个数据向量，简称数据；

网格：对地下待成像空间进行离散，离散后网格包含大量的格子，每个格子是一个多面体（如正六面体、不规则四面体等），每个格子包含的空间内有均一的电阻率；

模型：是对地下电阻率（或电导率）结构的描述，此处指狭义的模型，即网格中每个格子只有一个电阻率值，将所有格子的电阻率值取出按照一定顺序排列成一个向量，称为模型向量，简称模型；

观测模式：对地下空间施加天然源或人工源激励的方法，以及在空中、地面、地下测量电磁场的方法，具体包括：激发源的种类、参数选择和布设方法，接收器的种类、参数选择和布设方法以及野外观测的流程等；

正演：也称正向建模，指给定网络、模型和观测模式，通过数值模拟求解给出理论上应该得到的观测数据；

反演：也称反向建模，指给定观测模式、观测数据和网格，推算出地下空间电阻率模型的过程，与成像含义相同；

正演网格和正演模型：指计算正演响应数据时采用的求解数值解的网格和相应模型；

反演网格和反演模型：指用于容纳反演模型的网格和相应的反演成像结果模型。

接下来，对本申请实施例中提供的技术方案进行介绍。

实施例1

本申请提供的多尺度多方法地质体电阻率三维成像方法，涉及地球物理电磁法勘探领域，主要用于地球物理勘探（物探）数据的处理和解释；电磁勘探在野外采集数据之后，得到的是地球物理场的数据，这样的数据往往无法直接用于确定矿产、地下水和目标地质体的位置和埋深，因此，往往需要使用反演方法将地球物理场数据转换为地下空间的模型，以便后续的地质解释。

基于上述应用场景，如图1所示，本申请实施例提供了一种多尺度多方法地质体电阻率三维成像方法，包括以下步骤：

Step1：对地下空间使用任意剖分方法进行离散，获取反演网格，并在反演网格中填充均匀的电阻率值，构成初始反演电阻率模型，并作为当前反演电阻率模型；

这里需指出，任意剖分方法包括规则或不规则形体剖分，更为具体地包括：层状剖分、正交立体单元剖分或/和四面体剖分；如图2所示，反演网格根据地质条件从左到右使用规则网格、层状网格、非规则网格和半规则网格进行混合形式的离散化；

Step2：对不同类型的观测数据进行分组，以便于每组观测数据能够使用同样的数值求解方法进行正演；

更为具体地，假设观测手段包括直流电磁法、大地电磁法和瞬变电磁法；观测数据包含有n个不同类型的观测数据集；每个数据集分别包含有，…，/>个观测数据，所有观测数据共有/>个，对每一种观测数据集内的/>个观测数据值进行分组，分为/>组，所有数据分为/>个数据子集；

如图2所示，将平原地区且地层成层性较好的观测数据作为第一组观测数据，可以使用层状模型半解析解计算；将山地观测数据，且山体内有一个已知孔洞，则使用非规则网格剖分，使用有限元等方法求解的观测数据分为第二组观测数据；将另一平原地区，但地下存在局部不均匀体的观测数据作为第三组观测数据，使用基于规则网络的有限差分等方法计算；分组的目的是将具有相似性的数据放在一起，便于使用同样的数值求解方法进行正演，这里的相似性可以包括位置的相似性、地质结构的相似性和电磁感应尺度的相似性；

Step3：根据每组观测数据对应的观测手段、正演网格和当前反演电阻率模型，采用相应的数值求解方法进行正演，计算每组观测数据对应的理论正演响应数据，并将理论正演响应数据与实际观测数据进行求差，获取数据拟合差；

更为具体地，数据分组之后，每一组数据都会被给定一套对应的正演方法，正演方法包括：正演网格和正演模型、对应的观测模式和数值求解方法，共有套独立的计算方法；这/>套计算方法各自有对应的正演网格和正演模型；由于正演网格仅需要考虑某一组的数据，因此，每个正演网格仅需要考虑与其相关的空间覆盖范围、网格划分密度等等，更容易进行针对性优化；

这里需指出，初始反演电阻率模型和当前反演电阻率模型都定义在反演网格上，因此，为了在正演网格上计算当前反演电阻率模型的理论响应值，需要将当前反演电阻率模型映射到正演网格上；这一转换过程可以使用不同的方法，比如体积加权平均法，即为每个正演网格单元考虑与其有相交关系的反演网格单元，将相交体积作为权重，将反演网格单元内电阻率按照权重对正演网格单元内的平均电阻率进行贡献（如电阻率跨越多个数量级则使用对数电阻率）；得到正演网格和正演模型后，再配合该观测数据分组对应的观测手段，即可按照各分组各自的计算方法进行数值求解，输入正演方法后，计算出理论正演响应数据，与实际观测数据求差，得到数据拟合差；如果反演网格与正演网格不使用同一坐标系，在电阻率平均和正演结果汇总时还需要进行坐标转换；

Step4：若数据拟合差小于等于预设拟合差，则判定当前反演电阻率模型可以充分反映实测数据，反演步骤Step1终止，将当前反演电阻率模型用于地质体电阻率三维成像；若数据拟合差大于预设拟合差，则基于分组后的正演网格和正演模型计算m个缩减雅各比矩阵，转至步骤Step5；即缩减的雅各比矩阵的行数是分组内观测数据值的个数，列数是分组对应的正演网格的单元个数；

Step5：基于缩减雅各比矩阵和映射矩阵，获得完整雅各比矩阵的等效表示，完成完整雅各比矩阵与任意向量的乘积；

更为具体地，基于映射矩阵，获取正演模型与反演模型的电阻率对应关系，进而计算正演数据响应对当前反演电阻率模型内电阻率值的偏导数，获取各分组对应的完整雅各比矩阵；

Step6：使用完整雅各比矩阵和数据拟合差计算当前反演电阻率模型的更新量，更新当前反演电阻率模型，转至Step3，直至得到较好拟合实测数据的最终模型，或者达到最大允许迭代次数后反演过程强行终止；将最终的反演电阻率模型用于电阻率三维成像。

Step6具体包括以下步骤：

本申请中缩减雅各比矩阵都是在各个分组对应正演网格和正演模型上计算的，不能直接用于更新反演模型，所以需要再将正演网格上的缩减雅各比矩阵映射到反演网格上，结合数据拟合差信息，获取更新后的反演模型；

这里需指出，更新过程通过最优化目标函数实现，目标函数最优化问题通常使用牛顿、拟牛顿、高斯牛顿和共轭梯度等迭代方法求解；

在实际计算中，更新正演模型和更新反演模型的主要计算量最后都归结为雅各比矩阵（或其转置）与某个向量相乘；由于本申请在Step3中已经使用相交体积加权平均等方法计算了从反演网格向正演网格映射的矩阵，在Step6中将缩减雅各比矩阵映射到反演网格上获取反演网格上的完整雅各比矩阵，只需要将该映射矩阵转置之后作用于缩减雅各比矩阵，能得到等小的基于反演网格的完整雅各比矩阵；

实际上考虑完整雅各比矩阵与反演模型中的电阻率向量相乘的结果，这里并不需要把完整雅各比矩阵都算出来，只需要把映射矩阵作为雅各比矩阵-反演模型中的电阻率相乘中的一个步骤，即某个向量先与映射矩阵相乘，得到长度与正演模型单元个数相同的缩减向量，缩减向量再与缩减雅各比矩阵相乘，得到的结果等效于完整雅各比矩阵与某个向量相乘；由于缩减雅各比矩阵和映射矩阵分散存储于各个数据分组中，且互不相关，因此该步可以高效并行计算；

需指出，本申请中的映射矩阵反映的是正演网格单元与相邻反演网格单元之间的关系，本身是一个比较稀疏的矩阵；当正演网格单元显著大于相邻的反演网格单元时，映射矩阵可能出现局部稠密，增加内存占用；这时可以对映射矩阵进行行程编码、小波变换等压缩，近似处理以便于存储，但是计算时需要额外计算时间和空间解压缩；映射矩阵也可以通过单元中心点的三维插值获取，这一方法计算效率比较高，但是计算误差可能较大。

综上所述，本申请与现有技术相比，存在以下优势：

本申请对来自多个不同观测数据集进行分组，便于同组的观测数据使用相同的方式求解，而组和组之间可以不同，因此可以实现并行处理不同类型数据，使用不同数值求解方法处理不同尺度的观测数据。

本申请提出了分别为正演和反演设置反演网格和正演网格；反演网格和反演模型仅用于统筹不同观测数据以及反演结果，不用于数值计算，因此反演网格可以任意复杂（细网格、大范围和任意形状单元），不增加求解困难；具体数值求解由多个正演网格和正演模型完成，因此各个数据分组可以根据其具体条件和需求采用不同的最优方法，整体上可以达到性能的最优，不会陷入一个方法难以适应多种需求的困境。

本申请使用一个线性映射矩阵来处理从反演模型到正演模型的映射，以及从缩减雅各比矩阵到完整雅各比矩阵的映射；这两个映射互为转置关系，仅需要保存一个映射矩阵即可，得益于这一映射矩阵，反演过程中计算量最大最复杂的两个步骤（正演计算和雅各比矩阵与反演模型中的电阻率向量相乘）都可以在各个分组内完成，分组间互相独立便于并行求解，具有良好的可扩展性。

如果数据分组足够细，则正演网格单元数量非常小，因此有利于使用直接求解法计算麦克斯韦方程离散得到的线性方程组；对于小型线性方程组而言，直接求解比迭代解法计算速度快，计算精度高，内存消耗适中。

现有的电阻率模型的获取正反演共用的是一套网格的思路，本申请按照功能区分为反演网格和正演网格，且正演网格和正演模型可以根据数据分组为多个，互相独立；反演网格包含完整区域，有多个不同尺度、兼容各种地形地质和观测模式，每个正演网格仅负责一小部分数据、单一的尺度、单一的观测模式和某种特定的数值方法。

实施例2

如图3所示，本申请提供了一种多尺度多方法地质体电阻率三维成像系统，包括：

初始反演电阻率模型构建模块，用于对地下空间使用任意剖分方法进行离散后获取反演网格，并在反演网格中填充电阻率值，构建初始反演电阻率模型作为当前反演电阻率模型；

观测数据分组模块，用于对不同类型的观测数据进行分组，以便分组后的每组观测数据对应一个正演网格和正演模型，使用相同数值求解方法进行正演；

数据拟合差获取模块，用于将当前反演电阻率模型映射至正演网格，再根据每组观测数据对应的观测手段，采用相应的数值求解方法进行正演，计算每组观测数据对应的理论正演响应数据，并将理论正演响应数据与观测数据进行求差，获取数据拟合差；

反演电阻率模型判定模块，用于判断数据拟合差是否小于等于预设拟合差；若数据拟合差小于预设拟合差，则判定当前反演电阻率模型能够用于地质体电阻率三维成像；同时用于判断当前迭代次数是否达到最大迭代次数，若达到最大迭代次数则判定当前反演电阻率模型能够用于地质体电阻率三维成像；

反演模型的更新模块，用于若数据拟合差大于预设拟合差，则计算各正演模型对应的缩减雅各比矩阵；基于缩减雅各比矩阵和正演网格到反演网格之间的映射关系，完成定义在反演网格上的完整雅各比矩阵与任意向量的乘积操作，用以求取当前反演电阻率模型的总更新量向量；将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型；

其中，任意向量中包含有数据拟合差信息；缩减雅各比矩阵的行数是正演模型对应分组内观测数据的个数，列数为分组对应正演网格的单元个数。

进一步优选地，观测数据分组模块包括：基于观测手段的分组单元、基于相似性的分组单元和正演方法确定单元；

基于观测手段的分组单元用于根据获取观测数据采取的观测手段不同，将所述观测数据分为不同类型的观测数据集；

基于相似性的分组单元用于将每一组观测数据集内的观测数据根据位置的相似性、地质结构相似性和电磁感应尺度的相似性进一步划分；

正演方法确定单元用于根据划分后的每组观测数据确定对应的一个正演网格、正演模型以及一种数值求解方法。

进一步优选地，数据拟合差获取模块中将当前反演电阻率模型映射至正演网格的方法为：

将每个正演网格单元与反演网格单元相交的体积作为权重，构建映射矩阵；

将反演网格单元内电阻率按照映射矩阵对正演网格单元内的平均电阻率进行贡献。

进一步优选地，反演模型的更新模块包括：缩减雅各比矩阵求解单元、分组完整雅各比矩阵求解单元、更新向量求解单元和相加单元；

缩减雅各比矩阵求解单元用于求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

分组完整雅各比矩阵求解单元用于采用映射矩阵与各分组缩减雅各比矩阵相乘，获取各分组对应的完整雅各比矩阵；

更新向量求解单元用于综合各个分组完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更新向量；

相加单元用于将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型。

或

反演模型的更新模块包括：缩减雅各比矩阵求解单元、缩减向量求算单元、更新向量求解单元和相加单元；

缩减雅各比矩阵求解单元用于求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

缩减向量求算单元用于多次迭代搜索任意向量，并将任意向量与映射矩阵相乘，得到长度与正演模型单元个数相同的缩减向量；

更新向量求解单元用于采用缩减向量与各正演模型对应的缩减雅各比矩阵相乘，获得并综合各分组对应完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更细向量；

相加单元用于将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型。

在一些实施例中，当采集的观测数据来自平原地区且地层成层性超过预设成层性时，正演模型为层状模型，数值求解方法为解析解、半解析解求解方法；

当采集的观测数据来自山地且山体有孔洞时，反演网格为非规则网格剖分，数值求解方法为有限元求解方法；

当采集的观测数据来自平原地区且地下存在局部不均匀体时，反演网格为规则网格剖分，数值求解方法为有限差分、有限体积或积分方程方法。

应当理解的是，上述系统用于执行上述实施例中的方法，系统中相应的程序模块，其实现原理和技术效果与上述方法中的描述类似，该装置的工作过程可参考上述方法中的对应过程，此处不再赘述。

基于上述实施例中的方法，本申请实施例提供了一种电子设备。该设备可以包括：至少一个用于存储程序的存储器和至少一个用于执行存储器存储的程序的处理器。其中，当存储器存储的程序被执行时，处理器用于执行上述实施例中所描述的方法。

基于上述实施例中的方法，本申请实施例提供了一种计算机可读存储介质，计算机可读存储介质存储有计算机程序，当计算机程序在处理器上运行时，使得处理器执行上述实施例中的方法。

基于上述实施例中的方法，本申请实施例提供了一种计算机程序产品，当计算机程序产品在处理器上运行时，使得处理器执行上述实施例中的方法。

可以理解的是，本申请的实施例中的处理器可以是中央处理单元(centralprocessing unit，CPU)，还可以是其他通用处理器、数字信号处理器(digitalsignalprocessor，DSP)、专用集成电路(application specific integrated circuit，ASIC)、现场可编程门阵列(field programmable gate array，FPGA)或者其他可编程逻辑器件、晶体管逻辑器件，硬件部件或者其任意组合。通用处理器可以是微处理器，也可以是任何常规的处理器。

本申请的实施例中的方法步骤可以通过硬件的方式来实现，也可以由处理器执行软件指令的方式来实现。软件指令可以由相应的软件模块组成，软件模块可以被存放于随机存取存储器(random access memory，RAM)、闪存、只读存储器(read-only memory，ROM)、可编程只读存储器(programmable rom，PROM)、可擦除可编程只读存储器(erasable PROM，EPROM)、电可擦除可编程只读存储器(electrically EPROM，EEPROM)、寄存器、硬盘、移动硬盘、CD-ROM或者本领域熟知的任何其它形式的存储介质中。一种示例性的存储介质耦合至处理器，从而使处理器能够从该存储介质读取信息，且可向该存储介质写入信息。当然，存储介质也可以是处理器的组成部分。处理器和存储介质可以位于ASIC中。

在上述实施例中，可以全部或部分地通过软件、硬件、固件或者其任意组合来实现。当使用软件实现时，可以全部或部分地以计算机程序产品的形式实现。所述计算机程序产品包括一个或多个计算机指令。在计算机上加载和执行所述计算机程序指令时，全部或部分地产生按照本申请实施例所述的流程或功能。所述计算机可以是通用计算机、专用计算机、计算机网络、或者其他可编程装置。所述计算机指令可以存储在计算机可读存储介质中，或者通过所述计算机可读存储介质进行传输。所述计算机指令可以从一个网站站点、计算机、服务器或数据中心通过有线(例如同轴电缆、光纤、数字用户线(DSL))或无线(例如红外、无线、微波等)方式向另一个网站站点、计算机、服务器或数据中心进行传输。所述计算机可读存储介质可以是计算机能够存取的任何可用介质或者是包含一个或多个可用介质集成的服务器、数据中心等数据存储设备。所述可用介质可以是磁性介质，(例如，软盘、硬盘、磁带)、光介质(例如，DVD)、或者半导体介质(例如固态硬盘(solid state disk，SSD))等。

本领域的技术人员容易理解，以上所述仅为本申请的较佳实施例而已，并不用以限制本申请，凡在本申请的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本申请的保护范围之内。

Claims (10)

1.一种多尺度多方法地质体电阻率三维成像方法，其特征在于，包括以下步骤：

S1：对地下空间使用任意剖分方法进行离散后获取反演网格，并在反演网格中填充电阻率值，构建初始反演电阻率模型作为当前反演电阻率模型；

S2：对不同类型的观测数据进行分组，以便分组后的每组观测数据对应一个正演网格和正演模型，使用相同数值求解方法进行正演；

S3：将当前反演电阻率模型映射至正演网格，再根据每组观测数据对应的数值求解方法进行正演，计算每组观测数据对应的理论正演响应数据，并将理论正演响应数据与观测数据进行求差，获取数据拟合差；

S4：若数据拟合差小于等于预设拟合差，则判定当前反演电阻率模型能够用于地质体电阻率三维成像；否则转至S5；

S5：基于缩减雅各比矩阵和正演网格到反演网格之间的映射关系，完成定义在反演网格上的完整雅各比矩阵与任意向量的乘积操作，求取当前反演电阻率模型的总更新向量；其中，所述缩减雅各比矩阵的行数是正演模型对应分组内观测数据的个数，列数为分组对应正演网格的单元个数；所述任意向量中包含有数据拟合差信息；

S6：将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型，转至S3，直至数据拟合差小于等于预设拟合差或者迭代次数达到最大迭代次数，判定当前反演电阻率模型能够用于地质体电阻率三维成像。

2.根据权利要求1所述的多尺度多方法地质体电阻率三维成像方法，其特征在于，步骤S2具体包括以下步骤：

根据获取观测数据采取的不同观测手段，将所述观测数据分为不同类型的观测数据集；

再将每一组观测数据集内的观测数据根据位置的相似性、地质结构相似性和电磁感应尺度的相似性进一步划分；

根据划分后的每组观测数据确定对应的正演网格和正演模型以及不同的数值求解方法。

3.根据权利要求1或2所述的多尺度多方法地质体电阻率三维成像方法，其特征在于，将当前反演电阻率模型映射至正演网格的方法为：

将每个正演网格单元与反演网格单元相交的体积作为权重，构建映射矩阵；

将反演网格单元内电阻率按照映射矩阵对正演网格单元内的平均电阻率进行贡献。

4.根据权利要求3所述的多尺度多方法地质体电阻率三维成像方法，其特征在于，S5中求取当前反演电阻率模型的总更新量向量的方法为：

S5.1：求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

S5.2：采用映射矩阵与各分组缩减雅各比矩阵相乘，获取各分组对应的完整雅各比矩阵；

S5.3：多次迭代搜索并综合各个分组完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更新向量；

或

S5.1：求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

S5.2：多次迭代搜索任意向量，并将任意向量与映射矩阵相乘，得到长度与正演模型单元个数相同的缩减向量；

S5.3：采用缩减向量与各正演模型对应的缩减雅各比矩阵相乘，获得各分组对应完整雅各比矩阵与任意向量相乘的结果；

S5.4：综合各分组对应完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更新向量。

5.根据权利要求1或2所述的多尺度多方法地质体电阻率三维成像方法，其特征在于，当采集的观测数据来自平原地区且地层成层性超过预设成层性时，所述正演模型为层状模型，所述数值求解方法为解析解或半解析解求解方法；

当采集的观测数据来自山地且山体有异常体时，所述反演网格为非规则网格剖分，所述数值求解方法为有限元求解方法；

当采集的观测数据来自平原地区且地下存在局部不均匀体时，所述反演网格为规则网格剖分，所述数值求解方法为有限差分、有限体积或积分方程方法。

6.根据权利要求3所述的多尺度多方法地质体电阻率三维成像方法，其特征在于，当所述反演网格和所述正演网格不使用同一坐标系时，获取所述正演模型与当前反演电阻率模型的电阻率对应关系时进行坐标转换。

7.一种多尺度多方法地质体电阻率三维成像系统，其特征在于，包括：

初始反演电阻率模型构建模块，用于对地下空间使用任意剖分方法进行离散后获取反演网格，并在反演网格中填充电阻率值，构建初始反演电阻率模型作为当前反演电阻率模型；

观测数据分组模块，用于对不同类型的观测数据进行分组，以便分组后的每组观测数据对应一个正演网格和正演模型，使用相同数值求解方法进行正演；

数据拟合差获取模块，用于将当前反演电阻率模型映射至正演网格，再根据每组观测数据对应的观测手段，采用相应的数值求解方法进行正演，计算每组观测数据对应的理论正演响应数据，并将理论正演响应数据与观测数据进行求差，获取数据拟合差；

反演电阻率模型判定模块，用于判断数据拟合差是否小于等于预设拟合差；若数据拟合差小于预设拟合差，则判定当前反演电阻率模型能够用于地质体电阻率三维成像；同时用于判断当前迭代次数是否达到最大迭代次数，若达到最大迭代次数则判定当前反演电阻率模型能够用于地质体电阻率三维成像；

反演模型的更新模块，用于若数据拟合差大于预设拟合差，则计算各正演模型对应的缩减雅各比矩阵；基于缩减雅各比矩阵和正演网格到反演网格之间的映射关系，完成定义在反演网格上的完整雅各比矩阵与任意向量的乘积操作，用以求取当前反演电阻率模型的总更新量向量；将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型；

其中，所述任意向量中包含有数据拟合差信息；所述缩减雅各比矩阵的行数是正演模型对应分组内观测数据的个数，列数为分组对应正演网格的单元个数。

8.根据权利要求7所述的多尺度多方法地质体电阻率三维成像系统，其特征在于，所述观测数据分组模块包括：基于观测手段的分组单元、基于相似性的分组单元和正演方法确定单元；

基于观测手段的分组单元用于根据获取观测数据采取的观测手段不同，将所述观测数据分为不同类型的观测数据集；

基于相似性的分组单元用于将每一组观测数据集内的观测数据根据位置的相似性、地质结构相似性和电磁感应尺度的相似性进一步划分；

正演方法确定单元用于根据划分后的每组观测数据确定对应的一个正演网格、正演模型以及一种数值求解方法。

9.根据权利要求7或8所述的多尺度多方法地质体电阻率三维成像系统，其特征在于，数据拟合差获取模块中将当前反演电阻率模型映射至正演网格的方法为：

将每个正演网格单元与反演网格单元相交的体积作为权重，构建映射矩阵；

将反演网格单元内电阻率按照映射矩阵对正演网格单元内的平均电阻率进行贡献。

10.根据权利要求9所述的多尺度多方法地质体电阻率三维成像系统，其特征在于，所述反演模型的更新模块包括：缩减雅各比矩阵求解单元、分组完整雅各比矩阵求解单元、更新向量求解单元和相加单元；

缩减雅各比矩阵求解单元用于求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

分组完整雅各比矩阵求解单元用于采用映射矩阵与各分组缩减雅各比矩阵相乘，获取各分组对应的完整雅各比矩阵；

更新向量求解单元用于多次迭代搜索并综合各个分组完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更新向量；

相加单元用于将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型；

或

所述反演模型的更新模块包括：缩减雅各比矩阵求解单元、缩减向量求算单元、更新向量求解单元和相加单元；

所述缩减雅各比矩阵求解单元用于求取每组观测数据对应的理论正演响应数据对正演模型单元内电阻率的偏导数，获取缩减雅各比矩阵；

所述缩减向量求算单元用于多次迭代搜索任意向量，并将任意向量与映射矩阵相乘，得到长度与正演模型单元个数相同的缩减向量；

所述更新向量求解单元用于采用缩减向量与各正演模型对应的缩减雅各比矩阵相乘，获得并综合各分组对应完整雅各比矩阵与任意向量相乘的结果，获取当前反演电阻率模型的总更细向量；

所述相加单元用于将当前反演电阻率模型的总更新量向量与当前反演电阻率模型相加，获得更新后的反演电阻率模型。