技术领域
[0001] 本
申请涉及
地震波技术领域，具体而言，涉及一种基于垂直裂缝介质的地震波波场模拟方法及装置。
背景技术
[0002] 地震波正演模拟是模拟地震波在地球介质中的传播过程，以研究地震波传播特性与地下介质参数的关系。其中，数值模拟可以用于产生波传播的
可视化，从而便于直观的展示地下构造的真实形态。油气勘探实践表明，地质层裂缝是油气的储存空间或运移通道，地质层裂缝介质对地震波的影响的研究越来越受到关注。实际地质层中的裂缝形成受多种因素控制。由于上覆
地层的压
力使得
水平或低陡倾
角裂缝存在较少(即：水平裂缝)，大多数为高陡倾角裂缝(即：垂直裂缝)，因此，模拟垂直裂缝对地震波传播的影响更具意义。
[0003] 申请内容
[0004] 鉴于此，本申请
实施例的目的在于提供一种基于垂直裂缝介质的地震波波场模拟方法及装置，以模拟地质层的垂直裂缝对地震波传播的影响。
[0005] 第一方面，本申请实施例提供一种基于垂直裂缝介质的地震波波场模拟方法，所述方法包括：在已知地震波
信号、
震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的数值模拟表达式，其中，所述数值模拟表达式包括：
[0006]
[0007]
[0008] 以及
[0009]
[0010] 获取所述地质
密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述
采样时刻的值、所述时间采样率的值、所述空间
网格模型的大小，以及所述空间网格的长度值和宽度值；将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标
位置，以及所述空间网格的长度值和宽度值分别输入所述数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和
应力值。
[0011] 在上述实现过程中，在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，利用用于模拟所述垂直裂缝对所述地震波信号传播的影响的数值模拟表达式来获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值，继而能够根据所述速度值和所述应力值来模拟地质层的垂直裂缝对地震波传播的影响。
[0012] 基于第一方面，在一种可能的设计中，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的数值模拟表达式，包括：基于线性滑动理论，获取表征所述垂直裂缝位置处的边界条件；其中，所述边界条件为： 基于所述边界条件和虚拟网格，获取用于表征所述地震波信号在虚网格点处速度的第一表达式；获取表征应力和所述地震波信号在所述地质中的传播速度之间的关系的有限差分表达式；基于所述虚拟网格和所述旋转交错网格对所述有限差分表达式进行离散处理，得到离散化的差分表达式；基于所述第一表达式和所述离散化的差分表达式，获取所述数值模拟表达式。
[0013] 在上述实现过程中，由于对所述有限差分表达式进行离散处理的前提条件是连续，由于线性滑动理论可知，在垂直裂缝处速度是不连续的，因此，基于虚拟网格和旋转交错网格对所述有限差分表达式进行离散处理，得到离散化的差分表达式，并利用虚网格来对表征垂直裂缝位置处的边界条件处理，以得到表征所述地震波信号在虚网格点处速度的第一表达式，以显示的描述裂缝处的速度的不连续性，最后基于所述第一表达式和所述离散化的差分表达式得到所述数值模拟表达式，以解决上述问题。
[0014] 基于第一方面，在一种可能的设计中，基于所述虚拟网格和所述旋转交错网格对所述有限差分表达式进行离散处理，得到离散化的差分表达式，包括：利用所述虚拟网格和所述旋转交错网格对所述有限差分表达式中的速度和应力对所述第一差分方向和所述第二差分方向的偏导分别进行离散化处理，得到第一离散化的表达式；利用旋转交错网格对所述有限差分表达式中的速度和应力对时间的偏导分别进行离散化处理，得到第二离散化的表达式；利用所述第一离散化的表达式和所述第二离散化表达式对所述有限差分表达式进行离散处理，得到所述离散化的差分表达式。
[0015] 在上述实现过程中，为了方便求解所述有限差分表达式中的速度值和应力值，因此，通过上述方式对所述有限差分表达式进行离散处理，以方便求解所述有限差分表达式中的速度值和应力值。
[0016] 基于第一方面，在一种可能的设计中，获取表征应力和所述地震波信号在所述地质中的传播速度之间的关系的有限差分表达式，包括：获取表征所述应力和所述地震波信号在所述地质中的传播速度之间的关系的初始有限差分表达式；利用所述第一差分方向和所述第二差分方向对所述初始有限差分表达式中的差分方向进行处理，得到所述有限差分表达式。
[0017] 在上述实现过程中，为了方便计算地震波信号沿x轴和沿z轴的速度，以及沿x轴和z轴的应力，因此，利用旋转交错网格来改变所述初始有限差分表达式中的差分方向，继而方便计算地震波信号沿x轴和沿z轴的速度，以及沿x轴和z轴的应力。
[0018] 第二方面，本申请实施例提供一种基于垂直裂缝介质的地震波波场模拟方法，所述方法包括：在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的新的数值模拟表达式，其中，所述新的数值模拟表达式包括：
[0019]
[0020]
[0021]
[0022]
[0023]
[0024]
[0025] 以及
[0026]
[0027] 获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值；将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述新的数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0028] 在上述实现过程中，在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，利用用于模拟所述垂直裂缝对所述地震波信号传播的影响的新的数值模拟表达式来获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值，继而能够根据所述速度值和所述应力值来模拟地质层的垂直裂缝对地震波传播的影响。
[0029] 第三方面，本申请实施例提供一种基于垂直裂缝介质的地震波波场模拟装置，所述装置包括：第一模拟表达式获取单元，用于在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的数值模拟表达式；第一参数值获取单元，获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值；数值模拟单元，用于将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0030] 第四方面，本申请实施例提供一种基于垂直裂缝介质的地震波波场模拟装置，所述装置包括：第二模拟表达式获取单元，用于在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的新的数值模拟表达式，其中，所述新的数值模拟表达式包括；第二参数值获取单元，用于获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值；第二数值模拟单元，用于将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述新的数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0031] 第五方面，本申请实施例提供一种
电子设备，包括处理器以及与所述处理器连接的
存储器，所述存储器内存储
计算机程序，当所述计算机程序被所述处理器执行时，使得所述电子设备执行第一方面和第二方面所述的方法。
[0032] 第六方面，本申请实施例提供一种存储介质，所述存储介质中存储有计算机程序，当所述计算机程序在计算机上运行时，使得所述计算机执行第一方面和第二方面所述的方法。
[0033] 本申请的其他特征和优点将在随后的
说明书阐述，并且，部分地从说明书中变得显而易见，或者通过实施本申请实施例了解。本申请的目的和其他优点可通过在所写的说明书、
权利要求书、以及
附图中所特别指出的结构来实现和获得。
附图说明
[0034] 为了更清楚地说明本申请实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，应当理解，以下附图仅示出了本申请的某些实施例，因此不应被看作是对范围的限定，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。
[0035] 图1为本申请实施例提供的基于垂直裂缝介质的地震波波场模拟方法的流程示意图。
[0036] 图2为本申请实施例提供的地震波和垂直裂缝的场景图。
[0037] 图3为本申请实施例提高的旋转交错网格的网格示意图。
[0038] 图4为本申请实施例提供的虚拟网格示意图。
[0039] 图5为本申请实施例提供的表征虚网格和实网格重叠的示意图。
[0040] 图6为本申请实施例提供的另一种基于垂直裂缝介质的地震波波场模拟方法的
流程图。
[0041] 图7为本申请实施例提供的网格模型中的网格点的示意图。
[0042] 图8为本申请实施例提供的基于垂直裂缝介质的地震波波场模拟装置的结构示意图。
[0043] 图9为本申请实施例提供的另一种基于垂直裂缝介质的地震波波场模拟装置的结构示意图。
[0044] 图10为本申请实施例提供的地震波波场的一种传播示意图。
[0045] 图11为本申请实施例提供的地震波波场的另一种传播示意图。
[0046] 图12为本申请实施例提供的电子设备的结构示意图。
具体实施方式
[0047] 下面将结合本申请实施例中附图，对本申请实施例中的技术方案进行描述。
[0048] 应注意到：相似的标号和字母在下面的附图中表示类似项，因此，一旦某一项在一个附图中被定义，则在随后的附图中不需要对其进行进一步定义和解释。同时，在本申请的描述中，术语“第一”、“第二”等仅用于区分描述，而不能理解为指示或暗示相对重要性。
[0049] 请参照图1，图1为本申请实施例提供的一种基于垂直裂缝介质的地震波波场模拟方法，所述方法包括步骤：S100、S200和S300。
[0050] S100：在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的数值模拟表达式，其中，所述数值模拟表达式包括：
[0051]
[0052]
[0053]
[0054] 以及
[0055]
[0056] 其中，t表示采样时刻，(m，n)表示空间网格点在空间网格模型中的坐标位置，Δt表示所述时间采样率，Δx表示所述空间网格模型中的空间网格的长度，Δz表示所述空间网格的宽度，NT表示所述垂直裂缝的切向屈服度、NS表示所述垂直裂缝的法向屈服度，λ、μ表示拉梅系数，ρ表示地质密度， 表示所述地震波信号在采样时刻为t+1/2、空间网格点的坐标位置为(m，n)的沿x轴方向的速度， 表示所述地震波信号在采样时刻为t-1/2、空间网格点的坐标位置为(m，n)的沿z轴方向的速度， 表示所述地震波信号在采样时刻为t+1、空间网格点的坐标位置为(m+1/2，n+1/2)的沿x轴方向的应力，表示所述地震波信号在采样时刻为t、网格点的坐标位置为(m+1/2，n+1/2)的沿z轴方向的应力， 表示所述地震波信号在采样时刻为t、网格点的坐标位置为(m+1/2，n+1/2)的沿xz轴方向的应力。
[0057] S200：获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值；
[0058] S300：将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0059] 下面对上述方法进行详细阐述。
[0060] 在本实施例中，所述地震波信号为雷克子波，在其他实施例中所述地震波信号可以为高斯子波等信号波、瑞克子波等，请参考图2，在本实施例中，建立包含x轴和z轴的二维
坐标系，获取所述地震波信号开始传输的位置在所述二维坐标系中的坐标位置，即震源坐标(图中的点)，以及获取垂直裂缝(图中的线)在所述二维坐标系中的坐标位置，在本实施例中，所述震源坐标为(x＝750m，z＝750m)，所述垂直裂缝的坐标为(x＝1000m)，所述震源坐标和所述垂直裂缝的坐标根据实际需求设定，在此不做限制；在本实施例中，所述地震波信号的传播速度包括地震波的横
波速度和地震波的纵波速度，所述地震波的横波速度为1400m/s，地震波的纵波速度为2800m/s，在其他实施例中，所述地震波信号的传播速度根据实际需求设定，在此不做限制。
[0061] 作为一种实施方式，S100包括：获取表征所述应力和所述地震波信号在所述地质中的传播速度之间的关系的初始有限差分表达式；其中，所述初始有限差分表达式为：其中，vx表示在地质中传播的地震波信号沿x轴方向的速度分
量，vz表示表示在地质中传播的地震波信号沿z轴方向的速度分量，τxx表示沿x轴方向的应力分量，τxz表示沿xz轴方向的应力分量， 表示vx对时间t的偏导、 表示τxx对t的偏导， 表示应力τxz对x的偏导。由于获取所述初始有限差分表达式的方式为本领域技术人员所熟知的内容，因此在此不再赘述。
[0062] 为了求解所述初始有限差分表达式中的速度和应力的值，因此，在获取到所述初始有限差分表达式之后，利用第一差分方向和第二差分方向对所述初始有限差分表达式中的差分方向进行处理，获得表征应力和所述地震波信号在所述地质中的传播速度之间的关系的有限差分表达式；其中，所述有限差分表达式为：其中， 为旋
转交错网格的第一差分方向；为所述旋转交错网格的第二差分方向。
[0063] 在实际实施过程中，获得表征应力和所述地震波信号在所述地质中的传播速度之间的关系的有限差分表达式的步骤包括：A、B和C。
[0064] A：利用旋转交错网格定义所述第一差分方向和所述第二差分方向。请参照图3，下面对所述旋转交错网格的第一差分方向和所述旋转交错网格的差分方向进行介绍：
[0065] 本申请通过在预先建立的二维坐标系中建立网格模型，所述网格模型中有多个网格，定义所述第一差分方向 和所述第二差分方向 分别为沿着所述网格模型中的网格的两条对角线的上的方向，在所述网格模型大小相同的条件下，所述网格模型的长度和宽度的值越小，网格的数量越多，相反，网格的数量越少，网格的数量越多，因此，能够获取到的速度值和应力值越多，更有利于研究垂直裂缝对地震波的影响。
[0066] 作为一种实施方式，第一差分方向也可以为图3中的 方向，第二差分方向也可以为图3中的 方向。
[0067] B：在定义了所述第一差分方向和所述第二差分方向之后，利用所述第一差分方向和所述第二方向的线性组合来重新表示所述初始有限差分方程中的差分方向，得到新的差分方向表达式
[0068] C：利用所述新的差分方向表达式替换所述初始有限差分表示式中的差分方向，获得表征应力和所述地震波信号在所述地质中的传播速度之间的关系的有限差分表达式。
[0069] 为了求解所述有限差分表达式中的速度和应力值，因此，本申请需要利用差分算子代替微分算子，把计算区域离散化，以求取该区域每个网格
节点在不同时刻的速度值和应力值。
[0070] 下面对所述有限差分表达的离散过程进行详细阐述。
[0071] 请参照图4，图4中的垂直裂缝的左侧被设置为实网格，右侧被设置为虚网格，图4中的 为所述地震波信号在虚网格点坐标位置为(m+1，n)的速度；(v)m,n为所述地震波信号在实网格点坐标位置为(m，n)的速度。
[0072] 利用所述虚拟网格和所述旋转交错网格对所述有限差分表达式中的速度和应力对所述第一差分方向和所述第二差分方向的偏导分别进行离散化处理，得到第一离散化的表达式。
[0073] 其中，所述第一离散化的表达式包括：其中， 表示网格点位置为(m，n+1/2)处的速度(vi)m,n+1/2对 的偏导， 表
示网格点位置为(m，n+1/2)处的速度(vi)m,n+1/2对 的偏导， 表示网格点位置为(m，n+1/2)处的应力(τij)m,n+1/2对 的偏导，其中，vi包括vx和vz两种情况。
[0074] 利用旋转交错网格对所述有限差分表达式中的速度和应力对时间的偏导分别进行离散化处理，得到第二离散化的表达式。其中，所述第二离散化的表达式包括：其中，i和j均可以为x和z。例如，τij可以表示τxz也可以表示τxx
或τzz， 表示表示网格点位置为(m,n)处的应力(τij)m,n在第t时刻对时间的偏导。
[0075] 在获取到所述第一离散化表达式和所述第二离散化表达式之后，将所述第一离散化的表达式和所述第二离散化表达式代入所述有限差分表达式中以进行离散处理，得到离散化的差分表达式。
[0076] 其中，所述离散化的差分表达式包括：
[0077]
[0078]其中，Δt为时间采样率。其中，时间采用率Δt越小，模拟
精度越高，计算复杂的也越大。
[0079] 基于线性滑动理论，获取表征所述垂直裂缝位置处的边界条件；其中，所述边界条件为： 和 其中，“+”表示所述垂直裂缝的右侧，“-”表示所述垂直裂缝的左侧。
[0080] 其中，NT表示垂直裂缝的切向屈服度，Ns表示垂直裂缝的法向屈服度，vx+表示地震-波信号在垂直裂缝的左侧沿x轴方向的传播速度，vx表示地震波信号在垂直裂缝的右侧沿x轴方向的传播速度，vz+表示地震波信号在垂直裂缝的左侧沿z轴方向的传播速度，vz-表示地震波信号在垂直裂缝的右侧沿x轴方向的传播速度， 表示在垂直裂缝的左侧沿x轴方向的应力τxx+对时间t的偏导； 表示在垂直裂缝的右侧沿xz轴方向的应力τxz-对时间t的偏导。
[0081] 基于所述边界条件和虚拟网格，获取用于表征所述地震波信号在虚网格点处速度的第一表达式，其中，所述第一表达式包括：其中， 为所述地
震波信号在虚网格点坐标位置为(m，n)的沿x轴方向的速度。
[0082] 其中，获取用于表征所述地震波信号在虚网格点处速度的第一表达式的步骤可以按照如下方式实施：
[0083] 请参照图5，由于实网格和虚网格重叠在一起，因此，网格模型中的每个网格点的速度即可以用虚网格处的速度表示又可以用实网格处的速度表示，其中，图5中的裂缝的右侧用“+”表示，左侧用“-”表示，当需要计算裂缝处(即x方向)的速度时，从右侧看，速度从裂缝(图5中的波浪线)右侧传递到左侧时，由于地震波信号需要穿过裂缝，所以速度不连续，因此需要利用到虚网格处的速度，可以看作是从实网格到虚网格，反方向类似，继而根据中心差分准则，得到切向的速度差分表达式 根据算术平均法则，得到速度平均表达式 在需要计算裂缝左右两侧(即z方
向)的速度时，由于z方向没有裂缝，因此，地震波信号在z方向的速度为连续的，不需要利用虚网格，继而根据中心差分准则，得到法向的速度差分表达式
[0084] 由于所述初始有限差分表达式中的切向和法向的速度-应力表达式为：继而通过将所述切向的速度差分表达式、所述速度平均表达
式、法向的速度差分表达式以及所述速度-应力表达式，分别代入所述边界条件中经过简单的化简可得表征所述地震波信号在虚网格点处速度的第一表达式。
[0085] 在获取到了所述第一表达式和所述离散化的差分表达式之后，将所述第一表达式代入所述离散化的差分表达式中进行化简处理，获取所述数值模拟表达式。
[0086] S200：获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值。
[0087] 其中，所述地质密度的值ρ＝2600kg/m3，所述地质密度的值根据实际情况设定，所述拉梅系数的值λ＝0.2、μ＝0.8，拉梅系数的取值为经验值，在本实施例中，所述切向屈服度的值NT＝1.27*10-9，所述法向屈服度的值Ns＝2.69-9，NT和Ns的取值为经验值，在其他实施例中，NT和Ns的取值也可以为其他值，其中，所述空间网格模型的大小为1500m*1500m，所述网格模型中的每个空间网格的长度值为5m，宽度值为5m，其中，网格的长度值和宽度值不一定相等，在其他实施例中，所述网格模型的大小以及网格的长度值和宽度值也可以设置为其他值，例如，网格长度值为6m，宽度值为5m；在本实施例中，采样时刻t＝40ms，采样时刻也可以为20ms，30ms等，时间采用率Δt＝0.0001S，在其他实施例中，时间采样率的值根据实际需求设置，其中，时间采用率Δt越小，模拟精度越高，计算复杂的也越大。
[0088] S300：将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0089] 请参照图6，本申请实施例提供了另一种基于垂直裂缝介质的地震波波场模拟方法，所述方法包括：S400、S500和S600。
[0090] S400：在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的新的数值模拟表达式，其中，所述新的数值模拟表达式包括：
[0091]
[0092]
[0093]
[0094]
[0095]
[0096] 以及
[0097]
[0098]
[0099] 其中，t表示采样时刻，(m，n)表示空间网格点在空间网格模型中的坐标位置，Δt表示所述时间采样率，Δx表示所述空间网格模型中的空间网格的长度，Δz表示所述空间网格的宽度，NT表示所述垂直裂缝的切向屈服度、NS表示所述垂直裂缝的法向屈服度，λ、μ表示拉梅系数，ρ表示地质密度， 表示所述地震波信号在采样时刻为t+1/2、空间网格点的坐标位置为(m，n)的沿x轴方向的速度， 表示所述地震波信号在采样时刻为t-1/2、空间网格点的坐标位置为(m，n)的沿z轴方向的速度， 表示所述地震波信号在采样时刻为t+1、空间网格点的坐标位置为(m+1/2，n+1/2)的沿x轴方向的应力，表示所述地震波信号在采样时刻为t、网格点的坐标位置为(m+1/2，n+1/2)的沿z轴方向的应力， 表示所述地震波信号在采样时刻为t、网格点的坐标位置为(m+1/2，n+1/2)的沿xz轴方向的应力；C1和C2为有限差分系数。
[0100] S500：获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值。
[0101] S600：将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述新的数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0102] 利用有限差分进行数值模拟的原理是利用差分算子代替微分算子，把计算区域离散化，求取该区域每个网格节点在不同时刻的波场值。但是，当空间网格较大、震源主频过高时，会引起
波形的畸变，使波场发生较严重的数值频散。有限差分计算中数值频散压制思路主要为提高有限差分阶数，阶数越高数值模拟精度越高，但提高差分阶数会增加计算量。其主要原理是数值离散时，使用泰勒展开，根据泰勒展开，使用的阶数越高，精度越高。
[0103] 因此，为了抑制数值频散，以进一步提高数值模拟精度，因此，S400中获取所述新的数值模拟表达式的方式为：与S100不同的是，在获取到所述第一离散化表达式和所述第二离散化表达式之后，将所述第一离散化的表达式和所述第二离散化表达式带入所述有限差分表达式中以进行离散处理的过程中，在本申请实施例中，利用旋转交错网格的空间四阶来表示所述有限差分表达式中速度和应力对于第一差分方向和第二差分方向的偏导，请参照图7，在需要计算黑色方
块处的速度时，利用黑色圆点处的速度来计算黑色方块处的速度，图7中的网格点包括虚拟网格点和实网格点。
[0104] 在需要计算黑色方块处的速度时，若采用旋转交错网格的空间二阶的表示方式，则需要利用黑色方框附近的四个黑点的速度v(m，n)、v(m+1，n)、v(m，n+1)和v(m+1，n+1)，继而利用v(m+1，n+1)-v(m，n)代表 方向的二阶空间差分，另外两个相减代表 方向的二阶空间差分。若采用旋转交错网格的空间四阶的表示方式，则还需要利用到外面的四个黑点处的速度v(m+2，n+2)、v(m-1，n-1)、v(m+2，n-1)和v(m-1，n+2)，继而可以利用(v(m+1，n+1)-v(m，n))+(v(m+2，n+2)-v(m-1，n-1))代表 方向的四阶空间差分，利用上述方式对所述有限差分表达式进行离散化处理，继而得到新的离散化的差分表达式：
[0105]
[0106] 表1.几种有限差分系数
[0107]
[0108] 从表1可知，由于本申请是利用空间四阶来表示所述有限差分表达式中速度和应力对于第一差分方向和第二差分方向的偏导，因此，C1＝1.125，C2＝-0.0416666。
[0109] 通过将S100中获取到所述地震波信号在虚网格点处速度的第一表达式代入所述新的离散化的差分表达式，得到所述新的数值模拟表达式。其中，S500和S600的具体实施方式与S200和S300相同，因此在此不再赘述。
[0110] 请参照图8，本申请实施例提供一种基于垂直裂缝介质的地震波波场模拟装置，所述装置包括：第一模拟表达式获取单元410，用于在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的数值模拟表达式。第一参数值获取单元420，获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值。数值模拟单元430，用于将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0111] 作为一种实施方式，所述第一模拟表达式获取单元410，包括：边界条件获取单元，用于基于线性滑动理论，获取表征所述垂直裂缝位置处的边界条件；虚拟网格点处的速度获取单元，用于基于所述边界条件和虚拟网格，获取用于表征所述地震波信号在虚网格点处速度的第一表达式；差分表达式获取单元，用于获取表征应力和所述地震波信号在所述地质中的传播速度之间的关系的有限差分表达式；离散处理单元，用于基于所述虚拟网格和所述旋转交错网格对所述有限差分表达式进行离散处理，得到离散化的差分表达式；模拟表达式获取子单元，用于基于所述第一表达式和所述离散化的差分表达式，获取所述数值模拟表达式。
[0112] 作为一种实施方式，所述离散处理单元，还用于利用所述虚拟网格和所述旋转交错网格对所述有限差分表达式中的速度和应力对所述第一差分方向和所述第二差分方向的偏导分别进行离散化处理，得到第一离散化的表达式；利用旋转交错网格对所述有限差分表达式中的速度和应力对时间的偏导分别进行离散化处理，得到第二离散化的表达式；以及利用所述第一离散化的表达式和所述第二离散化表达式对所述有限差分表达式进行离散处理，得到所述离散化的差分表达式。
[0113] 作为一种实施方式，所述差分表达式获取单元，还用于获取表征所述应力和所述地震波信号在所述地质中的传播速度之间的关系的初始有限差分表达式；以及利用所述第一差分方向和所述第二差分方向对所述初始有限差分表达式中的差分方向进行处理，得到所述有限差分表达式。
[0114] 请参照图9，本申请实施例提供另一种基于垂直裂缝介质的地震波波场模拟装置，所述装置包括：第二模拟表达式获取单元510，用于在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的新的数值模拟表达式。第二参数值获取单元520，用于获取所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值；第二数值模拟单元530，用于将所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型中的每个空间网格点的坐标位置，以及所述空间网格的长度值和宽度值分别输入所述新的数值模拟表达式，得到用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值。
[0115] 本实施例对的各功能单元实现各自功能的过程，请参见上述图1-图7所示实施例中描述的内容，此处不再赘述。
[0116] 请参照图10和图11，图10为利用本申请所提供的方法在仿真实验中得到的地震波的垂直分量的传播示意图，图11为利用本申请所提供的方法在仿真实验中得到的地震波的水平分量的传播示意图，其中，仿真参数为：地震波采用主频为f＝30hz的瑞克子波；该地震波的纵波波速为2800m/s；该地震波的横波速度为1400m/s；地质密度ρ＝2600kg/m3。模型的空间网格参数Δx＝Δz＝5m，时间步长是Δt＝0.0001s；切向屈服度NT＝1.27×10-9m/Pa，法向屈服度Ns＝2.69×10-9m/Pa。
[0117] 请参照图12，图为本申请实施例提供的一种电子设备100的结构示意图，电子设备100可以包括：存储器102、处理器101和通信总线，通信总线用于实现这些组件的连接通信。
电子设备100可以为电脑、智能手机、
个人数字助理(personal digital assistant，PDA)等。
[0118] 存储器102用于存储所述地质密度的值、所述拉梅系数的值、所述切向屈服度的值、所述法向屈服度的值、所述采样时刻的值、所述时间采样率的值、所述空间网格模型的大小，以及所述空间网格的长度值和宽度值，所述数值模拟表达式，所述新的数值模拟表达式以及本申请实施例提供的基于垂直裂缝介质的地震波波场模拟方法及装置对应的计算指令，其中，存储器102可以是，但不限于，
随机存取存储器(Random Access Memory，RAM)，
只读存储器(Read Only Memory，ROM)等。
[0119] 其中，处理器101用于读取并运行存储于存储器102中的计算机指令时，执行本申请实施例所提供的基于垂直裂缝介质的地震波波场模拟方法的步骤。处理器101可能是一种集成
电路芯片，具有信号的处理能力。上述的处理器101可以是通用处理器，包括
中央处理器(Central Processing Unit，CPU)、网络处理器(Network Processor，NP)等。
[0120] 此外，本申请实施例还提供了一种存储介质，在该存储介质中存储有计算机程序，当所述计算机程序在计算机上运行时，使得所述计算机执行本申请任一项实施方式所提供的方法。
[0121] 综上所述，本申请各实施例提出的基于垂直裂缝介质的地震波波场模拟方法及装置，在已知地震波信号、震源坐标、所述地震波信号的传播速度、垂直裂缝坐标的条件下，利用用于模拟所述垂直裂缝对所述地震波信号传播的影响的数值模拟表达式来获取用于模拟所述垂直裂缝对所述地震波信号传播的影响的速度值和应力值，继而能够根据所述速度值和所述应力值来模拟地质层的垂直裂缝对地震波传播的影响。
[0122] 在本申请所提供的实施例中，应该理解到，所揭露的装置和方法，也可以通过其它的方式实现。以上所描述的装置实施例仅仅是示意性的，例如，附图中的流程图和
框图显示了根据本申请的多个实施例的装置和方法。在这点上，流程图或框图中的每个方框可以代表一个模块、程序段或代码的一部分，所述模块、程序段或代码的一部分包含一个或多个用于实现规定的逻辑功能的可执行指令。也应当注意，在有些作为替换的实现方式中，方框中所标注的功能也可以以不同于附图中所标注的顺序发生。例如，两个连续的方框实际上可以基本并行地执行，它们有时也可以按相反的顺序执行，这依所涉及的功能而定。也要注意的是，框图和/或流程图中的每个方框、以及框图和/或流程图中的方框的组合，可以用执行规定的功能或动作的专用的基于
硬件的装置来实现，或者可以用专用硬件与计算机指令的组合来实现。
[0123] 以上所述仅为本申请的优选实施例而已，并不用于限制本申请，对于本领域的技术人员来说，本申请可以有各种更改和变化。凡在本申请的精神和原则之内，所作的任何
修改、等同替换、改进等，均应包含在本申请的保护范围之内。
