基于先验梯度约束的四维电阻率跨孔CT监测成像方法与流程

本发明涉及一种基于先验梯度约束的四维电阻率跨孔ct监测成像方法。

背景技术

近年来，四维电阻率监测方法已经广泛的应用于环境监测，工程建设和水文地球物理研究等领域。跨孔勘察技术提供了比地面勘测上更高的分辨率和深度穿透能力。其测量传感器布置于钻孔中，信号不随深度衰减，有效的避免了工作现场中随机噪声的干扰。四维跨孔电阻率ct监测方法结合了跨孔勘察和四维电阻率层析成像技术的优势，是一种精细化监测的好手段。

目前，四维电法监测的研究取得了一定的效果，但仍存在以下关键问题，尚未解决：

其一，由于监测区域的地质环境复杂，被监测过程是动态变化的，采集的电阻数据十分复杂。使得四维反演方程的病态程度高，极易出现假异常，对目标体的定位不准确。

其二，现有的四维电阻率沿不同方向变化的敏感性是相同的。然而，在实际监测某些特定过程时，例如地下水运移，灌浆过程等，人们往往可以预先判断其运移的主要方向，这种先验信息并没有加到反演方程中。

技术实现要素：

本发明为了解决上述问题，提出了一种基于先验梯度约束的四维电阻率跨孔ct监测成像方法，本发明通过增加先验信息约束的方法改进反演方程，以突显电阻率在反演迭代过程中沿某一“主要方向”的时间演化特性，同时压制其他方向的演化特性，进而更有针对性的捕获目标体的时空位置、展布形态、运移路径和演化趋势。

为了实现上述目的，本发明采用如下技术方案：

一种基于先验梯度约束的四维电阻率跨孔ct监测成像方法，包括以下步骤：

利用多元地球物理手段进行探测，确定地质构造与环境状态信息，并根据探测结果划分重点监测区域；

根据场地条件在选定的重点监测区域布置钻孔，钻孔内布置电极，形成四维跨孔电阻率ct的观测空间；

在观测空间内确定正反演模型的范围，依据需要的精度和时间要求，确定电极间距以及正反演所建模型的网格大小；

开展全空间四维电阻率跨孔ct方法进行监测，构造四维电阻率反演目标函数与反演方程；

以探测的地质构造为指导，确定反演过程的初始模型，预先判断目标体的主要运移方向，作为先验信息，不断反演迭代得到一组电阻率反演结果，获得所监测区域动态变化的电导率图像。

进一步的，上述方法还包括根据获得的所监测区域动态变化的电导率图像，结合地质情况，解译目标体的时空位置、展布形态、运移路径和演化趋势。

进一步的，钻孔根据现场实际情况，选取合适地点，不需要严格的围成矩形，但尽量分布在目标体的四周。

进一步的，开展全空间四维电阻率跨孔ct精细监测时，利用不在同一直线上的三个或三个以上钻孔形成三维立体观测空间。

进一步的，在钻孔内布置测线，测线上均匀分布多个电极点，形成点源三维电场，对孔间地质结构进行多次覆盖采集和数据分析。

进一步的，所述四维电阻率反演目标函数为数据误差项、空间光滑约束项和时间光滑约束项之和。

更进一步的，所述数据误差项与实际观测数据与理论观测数据的差值向量相关。

更进一步的，所述空间光滑约束项与光滑约束矩阵以及模型参数在迭代中的模型参数增量向量的乘积相关。

更进一步的，在时间约束项中引入梯度项使得模型网格电阻率在反演迭代过程中沿三个方向具有不同的时间演化梯度特性。

进一步的，对模型进行独立的三维电阻率反演，由迭代一次后的结果，获得异常体的大致范围，对该范围以外的网格在三个方向上的梯度变化值都重设为0，令目标函数取极小，得到反演方程。

进一步的，根据地质信息以及地震探测结果，初步预测目标体运移的大致走向，作为方向梯度矩阵的权重代入反演方程，求解反演方程得到模型参数修正量，通过迭代得到反演结果。

与现有技术相比，本发明的有益效果为：

1.本发明提出了方向梯度矩阵的计算方法。通过反演迭代一次，确定探测对象内部异常目标体所在的空间区域，并对该区域加入方向梯度约束。目标体的运移扩散往往遵循一定的规律。方向梯度约束作为先验信息，可以从一定程度上对其动态过程进行表征，有利于更好的成像。

2.本发明提出了先验梯度信息权重的选择方式，由已知的地质信息和跨孔地震法探测的地下结构，对目标体的主要运移扩散方向作出判断，加大该方向的权重，使得反演结果中电阻率沿目标体主要运移方向的时间演化特性被描述地更为及时、更为准确。

3.本发明提出了基于时间先验梯度信息的四维电阻率反演方程，在时间约束项中加入梯度项，使得模型网格电阻率在反演迭代过程中沿三个方向具有不同的时间演化梯度特性，进一步加强了预测模型与参考模型间的关联。求解方程得到的一系列电导率时间图像，由此判断目标体的时空位置、展布形态、运移路径和演化趋势。

附图说明

构成本申请的一部分的说明书附图用来提供对本申请的进一步理解，本申请的示意性实施例及其说明用于解释本申请，并不构成对本申请的不当限定。

图1是基于先验梯度约束的四维电阻率跨孔ct监测成像方法的流程图。

图2是本发明进行数值模拟时使用的地电模型设计图。

图3是本发明进行数值模拟得到的反演结果。

具体实施方式：

下面结合附图与实施例对本发明作进一步说明。

应该指出，以下详细说明都是例示性的，旨在对本申请提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本申请所属技术领域的普通技术人员通常理解的相同含义。

需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本申请的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，当在本说明书中使用术语“包含”和/或“包括”时，其指明存在特征、步骤、操作、器件、组件和/或它们的组合。

在本发明中，术语如“上”、“下”、“左”、“右”、“前”、“后”、“竖直”、“水平”、“侧”、“底”等指示的方位或位置关系为基于附图所示的方位或位置关系，只是为了便于叙述本发明各部件或元件结构关系而确定的关系词，并非特指本发明中任一部件或元件，不能理解为对本发明的限制。

本发明中，术语如“固接”、“相连”、“连接”等应做广义理解，表示可以是固定连接，也可以是一体地连接或可拆卸连接；可以是直接相连，也可以通过中间媒介间接相连。对于本领域的相关科研或技术人员，可以根据具体情况确定上述术语在本发明中的具体含义，不能理解为对本发明的限制。

如图1所示，一种基于先验梯度约束的四维电阻率跨孔ct监测成像方法，包括以下步骤：

(1)在开展监测研究或试验之前，往往需要先通过多元地球物理手段(如地震波法、电法、电磁法等)进行普查，探明主要地质构造、地下水的径流状态等，在此基础上，根据探测结果划定和布置重点监测的区域。

(2)根据场地条件在选定的重点监测区域布置钻孔，钻孔内布置电极，形成四维跨孔电阻率ct的观测空间。

(3)依据步骤(2)中的观测空间，确定正反演模型的范围。依据实际所需的精度和时间要求，确定电极间距以及正反演所建模型的网格大小。

(4)开展全空间四维电阻率跨孔ct方法进行监测。构造四维电阻率反演目标函数及反演方程。

(5)以步骤(2)中探明的地质构造为指导，确定反演过程的初始模型，预先判断目标体的主要运移方向，作为先验信息代入携带(4)中方程，反演迭代得到一组电阻率反演结果。

(6)根据(5)中的成像结果，获得所监测区域动态变化的电导率图像，结合地质情况，解译目标体的时空位置、展布形态、运移路径和演化趋势。

步骤(2)中打的钻孔可以根据现场实际情况，选取合适地点，不需要严格的围成矩形，但尽量分布在目标体的四周。

步骤(4)中开展全空间四维电阻率跨孔ct精细监测，利用不在同一直线上的三个或三个以上钻孔形成三维立体观测空间，在孔内布置测线，测线上均匀分布多个电极点，形成点源三维电场，对孔间地质结构进行多次覆盖采集和数据分析。

步骤(4)中，所述全空间四维电阻率反演目标函数具体为：

φ(δm)＝||e^te||²+||wmδm||²+||wtct(m+δm)||²

等式右端第一项是数据误差项，e是实际观测数据与理论观测数据的差值向量。数据误差项中，定义采集的一系列数据集合d＝(d¹,d²,…,d^nc)^t。其中，dⁱ(i＝1,2,…,nc)表示第i次采集过程得到的电阻数据集，其中nc组完整数据集用于反演中。mⁱ(i＝1,2,…,nc)是一组模型参数，与dⁱ一一对应，m＝[m¹,...,m^nc]^t。e＝d-(g(m)+jδm)，d为实际观测数据，g代表正演过程。j为敏感度矩阵的扩展，j＝diag[j¹,j²,…,j^nc],j^k(k＝1,2,...,nc)为单个模型的敏感度矩阵，与m^k(k＝1,2,...,nc)一一对应。

右端第二项是空间光滑约束项，wm为光滑约束矩阵，δm是模型参数m在迭代中的模型参数增量向量。空间光滑约束使相邻网格电阻率差异极小，可表示如下：

式中，δmi为第i个网格的模型参数修正量，分别表示第i个网格的上、下、左、右、前、后的相邻网格的模型参数增量。将整个模型的光滑约束用矩阵形式表示，如下：

cmδm＝[r1,r2,...rner]^t

cm是对应于一个模型mⁱ(i＝1,2,…,nc)的空间约束矩阵，ner代表模型参数的个数，wm为cm在时间域上的扩展：

右端第三项即时间光滑约束项。其中wt为方向梯度光滑约束矩阵，ct为时间光滑约束矩阵，δm代表δm中的任一元素。时间光滑约束矩阵ct具体形式如下：

方向梯度光滑约束矩阵wt为对角块矩阵，wt＝diag(wt²,wt³,...,wt^nc)，wtⁱ(i＝2,3,…,nc)分别对应模型mⁱ(i＝2,3,…,nc)。

这里i是单位矩阵，η是单位矩阵的权重参数，是控制预测模型与参考模型偏差的光滑因子。梯度项其中gx，gy和gz分别代表模型mⁱ中x，y和z方向的方向梯度矩阵，表征模型电阻率沿x，y和z三个方向的变化梯度，α，β和γ为三个方向的权重系数。在时间约束项中引入梯度项使得模型网格电阻率在反演迭代过程中沿三个方向具有不同的时间演化梯度特性，通过调节α，β和γ的值，增大模型网格电阻率沿某一方向(如z方向)的梯度权重值，减小其他方向(如x、y方向)的梯度权重值。如令γ>α＝β，这将作为一种先验方向梯度约束信息使得反演迭代过程对网格电阻率沿z方向的时间演化特性更为敏感，最终使得时间推移反演成像结果中z方向的电阻率演化特性得到突显。

方向梯度矩阵计算方式如下：以gx为例，对任意网格，其电阻率值为r1，与沿x正方向的下一个网格，其电阻率值为r2，以(r2-r1)/r1作为该网格的梯度值，其中边界网格的梯度值设为0。将梯度值写入方向梯度矩阵中，长宽均为模型参数的个数。gy和gz同gx计算方式相同，不再累述。

需要说明的是，以上方式得到的矩阵还需要进一步处理。对模型mⁱ(i＝2,3,…,nc)进行独立的三维电阻率反演，由迭代一次后的结果，获得异常体的大致范围，对该范围以外的网格在三个方向上的梯度变化值都重设为0。

令以上目标函数取极小，得到如下反演方程：

所述步骤(5)中，地震法探测获得目标体的构造形态，，以此为基础，选取合适的初始模型用于(4)中反演方程。根据地质信息以及地震探测结果，初步预测目标体运移的大致走向，作为方向梯度矩阵的权重代入(4)中反演方程。求解反演方程得到δm，通过迭代得到反演结果。

a.地电模型图如附图2所示，反演区域设置为10m*10m*16m,背景电阻率为500ω·m。模型中添加了一个斜向下运移扩散的高阻体，其电阻率为1000ω·m。各阶段的异常体在模型中的位置和形态如附图2所示。

b.采集数据过程：正演使用矩形单元双线性插值的有限元方法，并采用自然边界条件。网格大小为0.5m*0.5m,电极间距设为1米,两根测线均布置16个电极。所用数据为bipole-bipole(am-bn)、pole–tripole(a-bmn和a-nmb)、dipole-dipole(ab-mn)三种传统的电极排列混合。

c.本算例中低阻体具有明显的垂向运动特征，换句话说，网格电阻率值在垂直方向上的变化将处于优势地位，因此，在时间加权矩阵中应增大对梯度分量gz的权重，本算例中时间加权矩阵表达式各系数的取值分别为μ＝0.5，η＝1.0，α＝β＝0.2，γ＝0.8。

d.经过反演处理之后的成像效果如附图3所示。该数值模拟表明，基于先验梯度约束的四维电阻率跨孔ct监测成像方法，能够有效降低反演方程的病态性，减少反演过程中的假异常，捕获目标体的时空位置、展布形态、运移路径和演化趋势。

以上所述仅为本申请的优选实施例而已，并不用于限制本申请，对于本领域的技术人员来说，本申请可以有各种更改和变化。凡在本申请的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本申请的保护范围之内。

上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。