CN109725345A - 一种初至波正演模拟方法及装置 - Google Patents

Abstract

本发明提供了一种初至波正演模拟方法及装置，其中，方法包括：确定深度域速度场，并对其进行离散化；根据深度域速度场设置震源子波信息以及观测系统信息；根据深度域速度场及观测系统信息，构造初至旅行时所用的程函方程，求取程函方程获取每个网格点上的初至时间场；按每个网格点上的初至时间场从小到大的顺序排序网格点序号；根据深度域速度场的最大速度值、最小的网格间距，计算得到波动方程数值正演时间间隔；根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新，实现初至波的正演模拟。本发明能够精确模拟初至波，极大的提高了初至波正演模拟的计算效率。

Description

一种初至波正演模拟方法及装置

技术领域

本发明涉及地球物理勘探技术，更具体地讲，涉及一种初至波正演模拟方法及装置。

背景技术

地震波场的数值模拟是研究地震波在复杂介质中传播规律的重要手段，同时也是一些反演问题(如全波形反演、逆时偏移等)的核心算法。地震波场正演模拟方法的研究不仅具有理论研究价值，同时还具有实际应用价值。

现有的地震波场正演模拟方法主要有：有限差分方法、有限体积方法、有限元方法等。由于有限差分方法相对较高的计算效率、简便的处理方式，应用最广；后两种方法，相对于有限差分方法对复杂地质构造具有更精确的描述能力，但由于计算效率较低，仅在特别复杂情况下应用。

现有的地震波场正演模拟方法均是基于波动方程实现的，在复杂近地表地区的研究中，通常只对初至波感兴趣，反射波等其他波型对近地表问题的研究贡献较小。由于复杂近地表区域存在起伏地表、复杂构造等复杂的地质构造，利用现有的数值算法进行正演模拟时往往需要特别小的空间网格及时间步长，从而会增加计算量，降低计算效率，另外，还会使初至波与反射波等波形混杂在一起，无法准确完整模拟初至波。现有的地震波场正演模拟方法极低的计算效率限制了初至波波形反演等近地表建模方法的应用。

因此，复杂近地表区域的研究需要构建一种计算效率更高的初至波正演模拟方法。

发明内容

本发明用于解决现有的地震波场正演模拟方法，在应用于近地表地区的研究时，存在计算效率低，初至波与反射波等波型混杂在一起的问题。

为了解决上述技术问题，本发明的第一方面提供一种初至波正演模拟方法，包括：

确定深度域速度场，并利用规则网格对深度域速度场进行离散化；

根据深度域速度场设置震源子波信息以及观测系统信息；

根据深度域速度场及观测系统信息，构造初至旅行时所用的程函方程，求取所述程函方程获取每个网格点上的初至时间场；

按每个网格点上的初至时间场从小到大的顺序排序网格点序号；

根据深度域速度场的最大速度值、最小的网格间距，计算得到波动方程数值正演时间间隔；

根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新，以实现初至波的正演模拟。

本发明进一步实施例中，深度域速度场表示为：

v(nx,ny,nz)，

其中，v为速度，nx、ny、nz分别为x、y、z三个方向离散后的网格点数，地表以上的网格点的速度用预定速度值填充。

本发明进一步实施例中，根据深度域速度场中的目标区地质体的尺寸设置震源子波信息。

本发明进一步实施例中，观测系统信息包括：炮点位置、检波点位置、炮点与检波点之间的关系以及最大记录时间。

本发明进一步实施例中，根据深度域速度场及观测系统构造的程函方程为：

程函方程满足如下边界条件：

t(x_s,y_s,z_s)＝0

其中，t(x,y,z)为炮点到空间位置(x,y,z)旅行时的离散初至时间场t_d(ix,iy,iz)的连续表示形式，v(x,y,z)为连续深度域速度场，(x_s,y_s,z_s)为炮点坐标位置；

采用有限差分方法求取程函方程，获取每个网格点上的初至时间场t_d(ix,iy,iz)，其中，ix、iy、iz为对应的x、y、z三个方向的网格点序号。

本发明进一步实施例中，计算得到波动方程数值正演时间间隔通过如下公式表示：

其中，α为小于1的系数，dtmax为最大时间间隔，dxyzmin为最小的网格间距，ndim为nx、ny、nz中不等于1的个数，vmax为深度域速度场的最大速度值，fdcoeff为包含有限差分系数的数组，nfdcoeff为所用的有限差分系数的个数。

本发明进一步实施例中，根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新的过程包括：

S610，对波动方程中的波场进行初始化；

S620，根据波场的时间迭代更新量及波动方程数值正演时间间隔，确定波场对应的正演时刻；

S630，根据正演时刻、初至时间场、排序后网格点序号以及震源子波信息，确定初至波空间范围；

S640，利用深度域速度场，基于波动方程对初至波空间范围的波场进行更新；

S650，获取检波点处的波场记录，同时时间迭代更新量加1；

S660，判断正演时刻是否到达最大记录时间，若未达到最大记录时间，则返回S620，否则，进入S670；

S670，将正演得到的初至波地震记录输出到文件。

本发明进一步实施例中，根据正演时刻、初至时间场、排序后网格点序号以及震源子波信息，确定初至波空间范围的过程包括：

将排序后网格点序号(ix,iy,iz)转换为一维索引序列order(i)；

按所述一维索引序列搜索得到一维索引序列序号idmin与idmax，使其满足：

其中，t_dl为一维索引形式的初至时间场，t_mod为正演时刻，τ为震源子波长度；

初至波空间范围Ω描述为：Ω＝{order[id]|idmin≤id≤idmax}。

本发明的第二方面提供一种初至波正演模拟装置，包括：确定模块，用于确定深度域速度场，并利用规则网格对深度域速度场进行离散化；

设置模块，用于根据深度域速度场设置震源子波信息以及观测系统信息；

第一计算模块，用于根据深度域速度场及观测系统信息，构造初至旅行时所用的程函方程，求取所述程函方程获取每个网格点上的初至时间场；

排序模块，用于按每个网格点上的初至时间场从小到大的顺序排序网格点序号；

第二计算模块，用于根据深度域速度场的最大速度值、最小的网格间距，计算得到波动方程数值正演时间间隔；

正演模块，用于根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新，以实现初至波的正演模拟。

本发明的第三方面提供一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，处理器执行所述计算机程序时实现前述任一初至波正演模拟方法实施例所述的方法。

本发明的第四方面提供一种计算机可读存储介质，计算机可读存储介质存储有执行计算机程序，计算机程序被处理器执行时实现前述任一初至波正演模拟方法实施例所述的方法。

本发明结合了基于射线理论的初至时间场与基于波动理论的波动方法正演模拟，通过初至时间场来约束每个时刻进行波场更新的范围，进而实现初至波的正演模拟，本发明适用于复杂近地表区域的研究，除了能够精确模拟初至波外，还能极大的提高了初至波正演模拟的计算效率。

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

附图说明

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

图1示出了本发明实施例初至波正演模拟方法的流程图；

图2示出了本发明实施例波场更新的流程图；

图3示出了本发明实施例初至波正演模拟装置的构成图；

图4示出了本发明实施例所用理论深度域速度场的示意图；

图5示出了本发明实施例正演模拟所得地震记录的示意图；

图6示出了现有技术地震波场正演模拟所得地震记录的示意图；

图7示出了图6地震记录与图5地震记录之差的示意图；

图8示出了本发明实施例初至波正演模拟方法在2.4s时的波场快照图；

图9示出了现有技术地震波场正演模拟方法在2.4s时的波场快照图；

图10示出了本发明初至波正演模拟方法与现有技术中地震波场正演模拟方法计算时间对比图。

具体实施方式

为了使本发明的技术特点及效果更加明显，下面结合附图对本发明的技术方案做进一步说明，本发明也可有其他不同的具体实例来加以说明或实施，任何本领域技术人员在权利要求范围内做的等同变换均属于本发明的保护范畴。

在本说明书的描述中，参考术语“一实施例”、“一具体实施例”、“一些实施例”、“例如”、或“一些实施方式”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不一定指的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。各实施例中涉及的步骤顺序用于示意性说明本发明的实施，其中的步骤顺序不作限定，可根据需要作适当调整。

如图1所示，图1示出了本发明实施例的初至波正演模拟方法的流程图，本实施例有效的结合了基于射线理论的初至时间场与基于波动理论的波动方法正演模拟，通过初至时间场来约束每个时刻进行波场更新的范围，进而实现初至波的正演模拟，本实施例适用于复杂近地表区域的研究，除了能够精确模拟初至波外，还能极大的提高了初至波正演模拟的计算效率。具体的，初至波正演模拟方法包括：

S100，确定深度域速度场，并利用规则网格对深度域速度场进行离散化；

S200，根据深度域速度场设置震源子波信息以及正演模拟所用的观测系统信息；

S300，根据深度域速度场信息及观测系统信息，构造初至旅行时所用的程函方程，求取所述程函方程获取每个网格点上的初至时间场；

S400，按每个网格点上的初至时间场从小大的顺序排序网格点序号；

S500，根据深度域速度场的最大速度值、最小的网格间距，计算得到波动方程数值正演时间间隔；

S600，根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新，以实现初至波的正演模拟。

本发明一些实施例中，上述步骤S100中，根据已有的地震、测井、地址等具有地质构造信息的资料，确定深度域速度场。深度域速度场使用三维数据进行存储，记为v(nx,ny,nz)，其中，nx、ny、nz分别为x、y、z三个方向离散后的网格点数，任一网格点的深度域速度场可以表示为v(ix,iy,iz)。地表以上的网格点的速度用预定速度值填充，预定速度值通常为0m/s。

本发明一些实施例中，上述步骤S200中，根据深度域速度场中的目标区地质体的尺寸设置震源子波信息，震源子波信息包含每个时间采样点上子波振幅、有效的子波长度τ、子波主频，其中，子波主频为15Hz到40Hz的雷克子波。

本发明一些实施例中，上述步骤S200中，根据深度域速度场的空间范围设定观测系统，观测系统具体包括：炮点位置、检波点位置、炮点与检波点之间的关系以及最大记录时间t_max。

本发明一些实施例中，上述步骤S300中，根据深度域速度场及观测系统构造的程函方程为

程函方程满足如下边界条件：

t(x_s,y_s,z_s)＝0 (式2)

其中，t(x,y,z)为炮点到空间位置(x,y,z)旅行时的离散初至时间场t_d(ix,iy,iz)的连续表示形式，v(x,y,z)为连续深度域速度场，(x_s,y_s,z_s)为炮点坐标位置。

采用有限差分方法求取程函方程，获取每个网格点上的初至时间场t_d(ix,iy,iz)，其中，ix、iy、iz为对应的x、y、z三个方向的网格点序号，地表以上的网格点的初至时间场设为极大值，一般设为1000s。

本发明一些实施例中，为了便于对符合空间范围的网格点进行标识，上述步骤S400中，采用快速排序方法按每个网格点上的初至时间场从小大的顺序排序网格点序号，将排序后的空间位置线性序号存储于一维数组order，数组大小为nx×ny×nz，数组中的第i个元素的值，用order(i)表示，其值的含义是空间位置的排序序号，与三维索引(ix,iy,iz)之间的关系为：

order(i)＝nz×(ix×ny+iy)+iz， (式3)

同样可将三维索引形式的初至时间场t_d(ix,iy,iz)转换为一维索引形式的时间场t_dl(ilinear)，两者满足：

t_dl(nz×(ix×ny+iy)+i_z)＝t_d(ix,iy,iz)， (式4)

排序后线性序号的一维数组order满足：

t_dl(order(i+1))≥t_dl(order(i))， (式5)

本发明一些实施例中，上述步骤S500中，波动方程数值正演的时间间隔dt为：

其中dtmax为最大时间间隔，dxyzmin为最小的网格间距，ndim为nx、ny、nz中不等于1的个数，vmax为深度域速度场的最大速度值，fdcoeff为包含有限差分系数的数组(该有限差分系数的数组可通过波动方程离散化后的频散关系与理论频散关系的最优拟合计算得到的)，nfdcoeff为所用的有限差分系数的个数。α为小于1的系数，一般选用0.9。

本发明一些实施例中，上述步骤S600中，所用的波动方程为一阶形式的常密度声波方程：

其中，p为压力波场，v为深度域速度场，v_x、v_y、v_z为质点位移速度场在x、y、z三个方向上的分量。当nx、ny、nz中某个维度为1时，将上式中对应的偏导设为零即可将上式退化为二维甚至是一维情况下的波动方程。f为震源项，具体的，震源项包括震源的空间位置、震源的时间域函数。

空间偏导可利用高阶交错网格有限差分进行数值计算，可利用算子表示为

其中u可以是p、v_x、v_y、v_z中的任意一种，算子D_x、D_y、D_z是有限长度的局部化算子，在其内部自动处理交错网格的不同情况。

本发明一些实施例中，如图2所示，上述步骤600根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行初至波的正演模拟的过程包括：

S610，对波动方程中的波场进行初始化。具体实施时，本步骤将波动方程中的压力波场p以及质点速度波场v_x、v_y、v_z在所有网格点上的值初始化为零，将时间迭代更新量it初始化为零。

S620，根据波场的时间迭代更新量it及波动方程数值正演时间间隔dt，确定波场对应的正演时刻t_mod＝it×dt。

S630，根据正演时刻t_mod、初至时间场、排序后网格点序号以及震源子波信息，确定初至波空间范围Ω。

按照排序后空间位置线性序号order搜索得到初至时间场一维索引序号idmin与idmax，使其满足：

其中，t_dl为一维索引形式的初至时间场，满足式(4)，t_mod为正演时刻，τ为震源子波长度；

初至波空间范围可用集合Ω描述为：

Ω＝{order[id]|idmin≤id≤idmax}， (式10)

S640，利用深度域速度场，基于波动方程对初至波空间范围Ω的波场进行更新。

具体实施时，波动方程中的时间偏导采用二阶精度的交错网格有限差分求解，其中压力波场p定义在规则时间采样点it上，质点速度场(v_x，v_y，v_z)定义时间采样中点it+1/2上。对于空间范围Ω，可应用(式3)将其转换回三维坐标索引表示的空间范围Ω_3D，对于(ix,iy,iz)∈Ω_3D的所有空间点，利用如下公式进行波场更新：

其中，p^it+1(ix,iy,iz)为时间采样点it+1、空间点(ix,iy,iz)的压力波场，p^it(ix,iy,iz)为时间采样点it、空间点(ix,iy,iz)的压力波场，v(ix,iy,iz)为空间点(ix,iy,iz)的深度域速度场，D_x、D_y、D_z为偏微分算子，分别为在时间采样点it+1/2、空间点(ix,iy,iz)、质点位移速度场在x、y、z三个方向上的分量，f^it+1/2为在时间采样点it+1/2的震源项， 分别为在时间采样点it+3/2、空间点(ix,iy,iz)、质点位移速度场在x、y、z三个方向上的分量，dt为波动方程数值正演时间间隔。

在更新完波场之后，通常需要将空间范围Ω_inner＝{order[id]|id＜idmin}内的波场进行清零处理。

S650，获取检波点处的波场记录。提取检波点对应位置处的压力波场，并用其更新波场记录，同时时间迭代更新量it加1。

S660，判断正演时刻t_mod是否到达最大记录时间t_max。若未达到最大记录时间，返回S620进行波场更新，否则，进入S670。

S670，将正演得到的初至波地震记录输出到文件。

基于同一发明构思，本发明实施例中还提供了一种初至波正演模拟装置，如下面的实施例所述。由于该装置解决问题的原理与初至波正演模拟方法相似，因此该装置的实施可以参见初至波正演模拟方法的实施，重复之处不再赘述。

该装置可以通过逻辑电路实现运行于智能终端，例如手机、平板电脑等设备中，或者以功能模块的方式由软件实现各部件的功能，运行于所述智能终端上，具体的，如图3所示，初至波正演模拟装置包括：

确定模块710，用于确定深度域速度场，并利用规则网格对深度域速度场进行离散化；

设置模块720，用于根据深度域速度场设置震源子波信息以及正演模拟所用的观测系统信息；

第一计算模块730，用于根据深度域速度场及观测系统信息，构造初至旅行时所用的程函方程，求取所述程函方程获取每个网格点上的初至时间场；

排序模块740，用于按每个网格点上的初至时间场从小到大的顺序排序网格点序号；

第二计算模块750，用于根据深度域速度场的最大速度值、最小的网格间距，计算得到波动方程数值正演时间间隔；

正演模块760，用于根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新，以实现初至波的正演模拟。

本实施例提供的初至波正演模拟装置，结合了基于射线理论的初至时间场与基于波动理论的波动方法正演模拟，通过初至时间场来约束每个时刻进行波场更新的范围，进而实现初至波的正演模拟，本实施例适用于复杂近地表区域的研究，除了能够模拟初至波外，还能极大的提高了初至波正演模拟的计算效率。

本发明一些实施例中，还提供一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，处理器执行所述计算机程序时实现前述任一初至波正演模拟方法实施例所述的步骤。初至波正演模拟方法的实施参见前述实施例，此处不再赘述。

本发明一些实施例中，还提供一种计算机可读存储介质，计算机可读存储介质存储有执行计算机程序，计算机程序被处理器执行时实现前述任一初至波正演模拟方法实施例所述的步骤。初至波正演模拟方法的实施参见前述实施例，此处不再赘述。

为了更清楚说明本发明技术方案，下面以一些具体实施例进行详细说明：

图4为测试所用的理论深度域速度场的示意图。该深度域速度场为2D速度场，深度方向最深7000m，采样间距为5m；横向共15000m，采样间隔为5m。图4中白色向下三角形为炮点位置，白色向上三角形为检波点位置。检波点横向起始位置为0m，横向间隔为40m，共260个检波点。基于图4所示深度域速度场利用本发明所述的正演模拟方法所得初至波记录如图5所示，利用现有技术地震波场正演模拟所得地震记录如图6所示，图7示出了图6地震记录与图5地震记录之差的结果。对比图4、图5及图6可知，本发明的初至波正演模拟方法能够精确地对初至波进行模拟。

如图8及图9所示，图8为利用本发明所述的初至波正演模拟方法在2.4s时的波场快照图，图9为利用现有技术地震波场正演模拟方法在2.4s时的波场快照图，通过对比图8及图9可以看出，本发明所述初至波正演模拟方法忽略了反射波等波形，仅模拟了初至波。

如图10所示，图10为本发明初至波正演模拟方法与现有技术中地震波场正演模拟方法计算时间对比图，通过对比可知，本发明所述的初至波正演模拟方法计算效率提高了12倍。

本发明结合了基于射线理论的初至时间场与基于波动理论的波动方法正演模拟，通过初至时间场来约束每个时刻进行波场更新的范围，进而实现初至波的正演模拟，本实施例除了能够模拟初至波外，还能极大的提高了初至波正演模拟的计算效率。

本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。

本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。

这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。

这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。

以上所述仅用于说明本发明的技术方案，任何本领域普通技术人员均可在不违背本发明的精神及范畴下，对上述实施例进行修饰与改变。因此，本发明的权利保护范围应视权利要求范围为准。

Claims (11)

1.一种初至波正演模拟方法，其特征在于，包括：

确定深度域速度场，并利用规则网格对深度域速度场进行离散化；

根据深度域速度场设置震源子波信息以及观测系统信息；

根据深度域速度场及观测系统信息，构造初至旅行时所用的程函方程，求取所述程函方程获取每个网格点上的初至时间场；

按每个网格点上的初至时间场从小到大的顺序排序网格点序号；

根据深度域速度场的最大速度值、最小的网格间距，计算得到波动方程数值正演时间间隔；

根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新，以实现初至波的正演模拟。

2.如权利要求1所述的方法，其特征在于，深度域速度场表示为：

v(nx,ny,nz)，

其中，v为速度，nx、ny、nz分别为x、y、z三个方向离散后的网格点数，地表以上的网格点的速度用预定速度值填充。

3.如权利要求1所述的方法，其特征在于，根据深度域速度场中的目标区地质体的尺寸设置震源子波信息。

4.如权利要求1所述的方法，其特征在于，观测系统信息包括：炮点位置、检波点位置、炮点与检波点之间的关系以及最大记录时间。

5.如权利要求1所述的方法，其特征在于，根据深度域速度场及观测系统构造的程函方程为：

程函方程满足如下边界条件：

t(x_s,y_s,z_s)＝0

其中，t(x,y,z)为炮点到空间位置(x,y,z)旅行时的离散初至时间场t_d(ix,iy,iz)的连续表示形式，v(x,y,z)为连续深度域速度场，(x_s,y_s,z_s)为炮点坐标位置；

采用有限差分方法求取程函方程，获取每个网格点上的初至时间场t_d(ix,iy,iz)，其中，ix、iy、iz为对应的x、y、z三个方向的网格点序号。

6.如权利要求1所述的方法，其特征在于，计算得到波动方程数值正演时间间隔通过如下公式表示：

其中，α为小于1的系数，dtmax为最大时间间隔，dxyzmin为最小的网格间距，ndim为nx、ny、nz中不等于1的个数，vmax为深度域速度场的最大速度值，fdcoeff为包含有限差分系数的数组，nfdcoeff为所用的有限差分系数的个数。

7.如权利要求1所述的方法，其特征在于，根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新的过程包括：

S610，对波动方程中的波场进行初始化；

S620，根据波场的时间迭代更新量及波动方程数值正演时间间隔，确定波场对应的正演时刻；

S630，根据正演时刻、初至时间场、排序后网格点序号以及震源子波信息，确定初至波空间范围；

S640，利用深度域速度场，基于波动方程对初至波空间范围的波场进行更新；

S650，获取检波点处的波场记录，同时时间迭代更新量加1；

S660，判断正演时刻是否到达最大记录时间，若未达到最大记录时间，则返回S620，否则，进入S670；

S670，将正演得到的初至波地震记录输出到文件。

8.如权利要求7所述的方法，其特征在于，根据正演时刻、初至时间场、排序后网格点序号以及震源子波信息，确定初至波空间范围的过程包括：

将排序后网格点序号(ix,iy,iz)转换为一维索引序列order(i)；

按所述一维索引序列搜索得到一维索引序列序号idmin与idmax，使其满足：

其中，t_dl为一维索引形式的初至时间场，t_mod为正演时刻，τ为震源子波长度；

初至波空间范围Ω描述为：Ω＝{order[id]|idmin≤id≤idmax}。

9.一种初至波正演模拟装置，其特征在于，包括：

确定模块，用于确定深度域速度场，并利用规则网格对深度域速度场进行离散化；

设置模块，用于根据深度域速度场设置震源子波信息以及观测系统信息；

第一计算模块，用于根据深度域速度场及观测系统信息，构造初至旅行时所用的程函方程，求取所述程函方程获取每个网格点上的初至时间场；

排序模块，用于按每个网格点上的初至时间场从小到大的顺序排序网格点序号；

第二计算模块，用于根据深度域速度场的最大速度值、最小的网格间距，计算得到波动方程数值正演时间间隔；

正演模块，用于根据深度域速度场、震源子波信息、观测系统信息、初至时间场、排序后网格点序号、波动方程数值正演时间间隔，基于波动方程进行波场更新，以实现初至波的正演模拟。

10.一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至8任一项所述的方法。

11.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有执行计算机程序，所述计算机程序被处理器执行时实现权利要求1至8任一项所述的方法。