加卸载响应比计算方法和系统与流程

本发明涉及地震预测技术领域，尤其是涉及一种加卸载响应比计算方法和系统。

背景技术：

近年来研究者发现，除了地震目录，使用地下水位、地壳形变等观测资料也能够发现到强震前异常变化，并推测可能与后续强震的发生有关。然而这些研究通常把潮汐因子作为判断加卸载的依据，因而难以达到规律性的认识。事实上，从微破裂产生、岩石物理力学性质改变到地震发生的过程中起决定作用的是有效剪切应力及其变化。

如何恰当地建立模型准确识别出微破裂的产生和扩展及对应的水位变化才是研究的关键问题所在。

技术实现要素：

有鉴于此，本发明的目的在于提供加卸载响应比计算方法和系统，通过简单方式计算加卸载响应比，提高计算结果的可靠性。

第一方面，本发明实施例提供了一种加卸载响应比计算方法,包括：

获取水位数据时间序列，所述水位数据时间序列包括按照多个预设时间点采集的水位数据；

对所述水位数据进行预处理和滤波操作，得到水位处理数据；

通过库伦应力准则、潮汐有效剪应力和构造有效剪应力，判断各个所述预设时间点的加卸载状态；

根据所述加卸载状态和所述水位处理数据，计算各个所述预设时间点的加卸载响应比。

结合第一方面，本发明实施例提供了第一方面的第一种可能的实施方式，其中，所述通过库伦应力准则、潮汐有效剪应力和构造有效剪应力，判断各个所述预设时间点的加卸载状态包括：

根据库伦应力准则，计算各个所述预设时间点对应的所述潮汐有效剪应力在所述构造有效剪应力上的投影值；

对各个所述投影值进行微分计算，得到变化率；

根据所述变化率的正负，判断各个所述预设时间点的加卸载状态。

结合第一方面，本发明实施例提供了第一方面的第二种可能的实施方式，其中，所述根据所述加卸载状态和所述水位处理数据，计算各个所述预设时间点的加卸载响应比包括：

对所述预设时间点为加载状态的水位处理数据进行求平均值计算，得到加载平均值；

对所述预设时间点为卸载状态的水位处理数据进行求平均值计算，得到卸载平均值；

将所述加载平均值与所述卸载平均值进行作商运算，得到所述预设时间点的加卸载响应比。

结合第一方面，本发明实施例提供了第一方面的第三种可能的实施方式，其中，在所述通过库伦应力准则、潮汐有效剪应力和构造有效剪应力，判断各个所述预设时间点的加卸载状态的步骤之前，还包括：

通过固体潮获取各个所述预设时间点对应的潮汐有效剪应力。

结合第一方面，本发明实施例提供了第一方面的第四种可能的实施方式，其中，所述对所述水位数据进行预处理和滤波操作，得到水位处理数据包括：

去除突跳点的水位数据；

对不完整的水位数据进行插值处理；

经带通滤波操作，得到水位处理数据。

结合第一方面，本发明实施例提供了第一方面的第五种可能的实施方式，其中，所述去除突跳点的水位数据包括：

判断当前所述预设时间点是否为突跳点；

当所述预设时间点为突跳点时，将所述预设时间点对应的水位数据进行删除。

结合第一方面，本发明实施例提供了第一方面的第六种可能的实施方式，其中，所述判断当前所述预设时间点是否为突跳点包括，重复执行以下处理，直至所述水位数据时间序列中的各个所述预设时间点对应的所述水位数据都被遍历：

计算当前预设时间点对应的水位数据与前一预设时间点对应的水位数据的差值；

若所述差值的绝对值超过预设范围，则所述当前预设时间点为突跳点，将所述突跳点对应的所述水位数据进行删除。

结合第一方面，本发明实施例提供了第一方面的第七种可能的实施方式，其中，所述对不完整的水位数据进行插值处理包括：

识别出所述水位数据时间序列中的所述不完整的水位数据对应的空缺时间点；

根据所述空缺时间点相邻所述预设时间点对应的所述水位数据的平均值，估算所述空缺时间点对应的所述不完整水位数据；

通过线性插值技术将所述不完整水位数据插入所述水位数据时间序列中，得到所述水位处理数据。

结合第一方面，本发明实施例提供了第一方面的第八种可能的实施方式，其中，所述根据所述变化率的正负，判断各个所述预设时间点的加卸载状态包括：

当所述变化率为正时，所述预设时间点为加载状态；

当所述变化率为负时，所述预设时间点为卸载状态。

第二方面，本发明实施例还提供一种加卸载响应比计算系统，包括：

数据获取单元，获取水位数据时间序列，所述水位数据时间序列包括按照多个预设时间点采集的水位数据；

数据处理单元，对所述水位数据进行预处理和滤波操作，得到水位处理数据；

状态判断单元，通过库伦应力准则、潮汐有效剪应力和构造有效剪应力，判断各个所述预设时间点的加卸载状态；

计算单元，根据所述加卸载状态和所述水位处理数据，计算各个所述预设时间点的加卸载响应比。

本发明实施例提供了一种加卸载响应比计算方法和系统，通过对获取的水位数据进行预处理、滤波操作，提高水位数据的质量，再根据固体潮变化获取潮汐有效剪应力，根据库伦应力准则、潮汐有效剪应力判断当前预设时间点的加卸载状态，分别计算出加载状态和卸载状态的水位平均值，将加载状态的水位平均值与卸载状态的水位平均值作商运算，得到加卸载响应比，通过简单有效的计算方式，得到更为可靠的计算结果，便于地震预测系统根据本发明实施例提供可靠的加卸载响应比得到更加准确的预测结果；

本发明的其他特征和优点将在随后的说明书中阐述，并且，部分地从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点在说明书、权利要求书以及附图中所特别指出的结构来实现和获得。

为使本发明的上述目的、特征和优点能更明显易懂，下文特举较佳实施例，并配合所附附图，作详细说明如下。

附图说明

为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例提供的加卸载响应比计算方法流程图；

图2为本发明实施例提供的发震断层面、滑移矢量以及构造和潮汐有效剪应力的叠加示意图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

目前,震前扩容是建立在块体岩石材料的本构关系之上,大量岩石实验和地震研究都能够为扩容模型提供支持。研究发现一旦观测确认了扩容过程，许多其他类型的前兆现象也可以利用这一扩容模型予以解释，然而地震波速(比)的测量精度却制约了这一模型的使用和发展。事实上，岩体的扩容和孔隙内流体的流动直接相关，其变化尺度足以使水位产生可测得的变化。这可能在一定程度上弥补了观测资料上的不足，从另一个角度延续和发展了扩容-扩散模型。

基于此，本发明实施例提供的一种加卸载响应比计算方法和系统，可以通过简单方式计算加卸载响应比，提高计算结果的可靠性,为扩容模型的广泛应用发展提供便利。

为便于对本实施例进行理解，首先对本发明实施例所公开的一种加卸载响应比计算方法进行详细介绍；

图1为本发明实施例提供的加卸载响应比计算方法流程图；

参照图1,加卸载响应比计算方法包括以下步骤:

步骤s110，获取水位数据时间序列，水位数据时间序列包括按照多个预设时间点采集的水位数据；

步骤s120，对水位数据进行预处理和滤波操作，得到水位处理数据；

步骤s130，通过库伦应力准则、潮汐有效剪应力和构造有效剪应力，判断各个预设时间点的加卸载状态；

步骤s140，根据加卸载状态和水位处理数据，计算各个预设时间点的加卸载响应比。

在实际应用的优选实施例中，将库仑应力变化结合震前扩容的概念判断水位的加卸载过程，得到加卸载响应比，通过对比加卸载过程的水位平均值变化，分析地震的迫近程度，计算方式简单，同时计算得到的响应比结果更加可靠；

进一步的，为了提高水位数据的质量，对水位数据进行预处理操作，包括去突跳和数据补全，上述实施例中的步骤s120包括：

步骤s210，去除突跳点的水位数据；

需要说明的是，裂隙是构造应力积累到一定程度后在岩石体内产生的，由裂隙引起的体积变化不超过正常弹性体积变形的0.2-2倍，因此发生瞬时的大幅突跳水位数据的不在裂隙引起的体积变化范围内，应予以去除，以提高响应比计算的可靠性；去除过程中，应首先判断当前的预设时间点是否为突跳点；当预设时间点为突跳点时，将预设时间点对应的水位数据进行删除，以提高水位数据的质量和稳定性；

其中，在一些可能的实施例中，可通过以下步骤，具体判断当前的预设时间点是否为突跳点，即重复执行以下处理，直至水位数据时间序列中的各个预设时间点对应的水位数据都被遍历：

计算当前预设时间点对应的水位数据与前一预设时间点对应的水位数据的差值；

若差值的绝对值超过预设范围，则当前预设时间点为突跳点，将突跳点对应的水位数据进行删除。

步骤s220，对不完整的水位数据进行插值处理；

这里，首先识别出水位数据时间序列中的不完整的水位数据对应的空缺时间点；再根据空缺时间点相邻预设时间点对应的水位数据的平均值，估算空缺时间点对应的不完整水位数据，并记录不完整水位数据；再通过线性插值技术将不完整水位数据插入水位数据时间序列中，得到水位处理数据；

步骤s230，经带通滤波操作，得到水位处理数据。

其中，步骤s230使用巴特沃斯(butterworth)滤波器对水位数据进行带通滤波处理，只保留周期为1-24h频段的水位数据，以去除非构造因素的干扰信号，提高水位数据的稳定性。

进一步的,在步骤s130之前，还包括：通过固体潮在特定滑移方向上引起的库仑应力变化，获取各个预设时间点对应的潮汐有效剪应力。

进一步的,上述实施例中的步骤s130还可用以下步骤实现：

步骤s310，根据库伦应力准则，计算各个预设时间点对应的潮汐有效剪应力在构造有效剪应力上的投影值；

步骤s320，对各个投影值进行微分计算，得到变化率；

步骤s330，根据变化率的正负，判断各个预设时间点的加卸载状态。

这里，当变化率为正时，预设时间点为加载状态；当变化率为负时，预设时间点为卸载状态。

在实际应用的过程中，根据岩石受载的应力-应变关系，当载荷积累到一定程度以后，产生的现象会与线弹性阶段不同，由于偏应力增加，岩石体内散布的裂隙会增加，与之相伴的是地下水位可能发生变化，由于变形的不可逆性，在循环载荷下，加卸载过程的水位就会出现差异。由于裂隙的产生反应了岩石本构关系的动态变化，因此这种差异可能与后续强震的发生密切相关；

需要说明的是本发明实施例通过判别岩石的体积动态变化，衡量水位变化，进而可通过上述变化情况对地壳活动异常进行监测；

根据库伦应力准则可知，由于裂隙是剪应力所致，因此我们采用库仑应力准则判断加卸载过程，其中，库仑破裂应力cfs定义为：

cfs＝τ-fσ(1)

公式(1)中f代表内摩擦系数，τ表示剪应力和σ表示法向应力，其中，公式(1)中的τ-fσ也被称为有效剪应力τe(zobacketal.,1987)，通常可以表示为构造有效剪应力和潮汐有效剪应力之和：

相对于构造应力的加载，如在地下10公里处，潮汐应力与构造应力相比低了约2-3个数量级，因此加载方向与构造应力方向相一致。如果假定观测点处的剪切运动方向与对应地震断层面上的滑移矢量方向相一致，具体如图2所示，构造有效剪应力和潮汐有效剪应力可以分别表示为：

图2(a)表示发震断层面和滑移矢量，构造和潮汐有效剪应力的叠加示意分别参考图2(b)和图2(c)。

另一方面，与潮汐应力的时间尺度(1天)相比，构造应力的变化却非常缓慢，二者相差超过4个数量级。因此，构造应力的变化率远小于潮汐应力的变化率，即：换句话说，在潮汐应力变化的时间尺度，构造应力几乎可看作不变的，这样，再由公式(2)可以将有效剪应力变化近似的由潮汐有效应力的变化表示，如下公式(4)所示：

其中，公式(4)中δτe可理解为潮汐有效剪应力在构造有效剪应力矢量方向上的投影值；

结合图2(b)可以看到，当构造有效剪应力和潮汐有效剪应力的夹角θ小于90°时，构造有效剪应力和潮汐有效剪应力互相加强，当构造有效剪应力和潮汐有效剪应力的夹角θ大于90°时，二者互相减弱。

进一步对上述投影值进行取微分计算，可以得到变化率g，如公式(5)所示：

其中，当g>0时为加载状态，当g<0时为卸载状态。

这里要指出的是，潮汐应力可以通过计算地壳弹性变形的方法获得，地壳变形可以用6个一阶微分方程表示，使用数值方法(runge-kutta)可以计算任意截面上的潮汐应力分量，再通过应力张量变换可以得到滑移面上的剪应力和正应力。

进一步的,根据岩石本构关系，当载荷积累到一定程度，偏离线弹性阶段以后，由于库仑破裂应力增加，岩石体内的微破裂会增加，从而导致体积变化，与之相伴的是地下水位发生改变。由于过程的不可逆性，按照凯塞(kaiser)效应，循环载荷下，加载和卸载过程的水位响应会出现差异。这种水位差异所反映的岩石介质本构关系的动态变化，可通过加卸载响应比数值体现，该响应与后续地震的发生具有关联性。

这里，根据上述步骤s140计算加卸载响应比，包括：

步骤s410，对预设时间点为加载状态的水位处理数据进行求平均值计算，得到加载平均值；

步骤s420，对预设时间点为卸载状态的水位处理数据进行求平均值计算，得到卸载平均值；

步骤s430，将加载平均值与卸载平均值进行作商运算，得到预设时间点的加卸载响应比。

在实际计算中，我们取一段时间内的平均水位作为响应量，即选取一定长度时间窗的水位数据时间序列，根据时间序列中的水位数据平均值作为响应量，以此研究加卸载过程中的水位变化，从而将加卸载响应比定义为公式(6)所示：

其中，n+和n-分别为加载状态和卸载状态的水位数据的记录总数，hi表示第i个水位数据值。时间窗内通常包含多个加、卸载过程，其时间窗长度与地震震级有关，一般震级越高，时间窗长度越长。

本发明实施例通过对获取的水位数据进行预处理、滤波操作，提高水位数据的质量，再根据固体潮变化获取潮汐有效剪应力，根据库伦应力准则、潮汐有效剪应力判断当前预设时间点的加卸载状态，分别计算出加载状态和卸载状态的水位平均值，将加载状态的水位平均值与卸载状态的水位平均值作商运算，得到加卸载响应比，通过简单有效的计算方式，得到更为可靠的计算结果，便于地震预测系统根据本发明实施例提供可靠的加卸载响应比得到更加准确的预测结果；

进一步的,本发明实施例还提供一种加卸载响应比计算系统，包括：

数据获取单元，获取水位数据时间序列，水位数据时间序列包括按照多个预设时间点采集的水位数据；

数据处理单元，对水位数据进行预处理和滤波操作，得到水位处理数据；

状态判断单元，通过库伦应力准则、潮汐有效剪应力和构造有效剪应力，判断各个预设时间点的加卸载状态；

计算单元，根据加卸载状态和水位处理数据，计算各个预设时间点的加卸载响应比。

本发明实施例提供的加卸载响应比计算系统，与上述实施例提供的加卸载响应比计算方法具有相同的技术特征，所以也能解决相同的技术问题，达到相同的技术效果。

本发明实施例所提供的加卸载响应比计算方法和系统的计算机程序产品，包括存储了程序代码的计算机可读存储介质，所述程序代码包括的指令可用于执行前面方法实施例中所述的方法，具体实现可参见方法实施例，在此不再赘述。

所属领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的系统和装置的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。

另外，在本发明实施例的描述中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。

所述功能如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：u盘、移动硬盘、只读存储器(rom，read-onlymemory)、随机存取存储器(ram，randomaccessmemory)、磁碟或者光盘等各种可以存储程序代码的介质。

在本发明的描述中，需要说明的是，术语“中心”、“上”、“下”、“左”、“右”、“竖直”、“水平”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。此外，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性。

本发明实施例还提供一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，处理器执行计算机程序时实现上述实施例提供的加卸载响应比计算方法的步骤。

本发明实施例还提供一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序，计算机程序被处理器运行时执行上述实施例的加卸载响应比计算方法的步骤。

最后应说明的是：以上所述实施例，仅为本发明的具体实施方式，用以说明本发明的技术方案，而非对其限制，本发明的保护范围并不局限于此，尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，其依然可以对前述实施例所记载的技术方案进行修改或可轻易想到变化，或者对其中部分技术特征进行等同替换；而这些修改、变化或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应所述以权利要求的保护范围为准。