基于预计算的煤矿井下二维波速分布反演方法

本发明涉及一种煤矿井下二维波速分布反演方法，属于矿震波速反演和智能算法。

背景技术：

1、目前，矿震波速反演或矿震波速成像是一种获得煤矿井下应力分布的有效方法，并已在井下冲击矿压预测上取得良好效果。由于检波器布设位置受限，目前主流的反演思路是将矿震波速反演问题转化为最优化问题求解，这一思路已取得很多的成果，但同时也带来应选择何种最优化方法、应该如何设计目标函数等问题。因此，有必要设计一种具备更快的反演速度和更高的可靠性的反演方法及系统。

技术实现思路

1、为了克服上述现有技术的不足之处，本发明提供一种基于预计算的煤矿井下二维波速分布反演方法，通过检索模型库，给出了其中符合目标函数的所有可能结果，具备更快的反演速度和更高的可靠性。

2、本发明是通过如下技术方案实现的：一种基于预计算的煤矿井下二维波速分布反演方法，包括二维波速分布模型库构建和反演两个过程：

3、其中所述二维波速分布模型库构建过程包括如下步骤：

4、步骤1：

5、步骤1-1：设定目标区域空间分辨率m、n：即待反演区域网格数，此区域以正方形剖分成m×n的网格，剖分所用正方形边长为dx＝dy。

6、设定空域波速最大值vgmax、最小值vgmin：目标区域以正方形进行剖分后，认为每个网格内的波速均匀分布，vgmax对应所有网格的波速可能达到的最大值，vgmin对应所有网格的波速可能达到的最小值，vgmax与vgmin的单位均为米/秒。

7、设定空域平均波速最大值最小值目标区域以正方形进行剖分后，认为每个网格内的波速均匀分布，空域平均波速为目标区域所有网格波速的平均值；定义目标区域的某个可能的波速分布对应波速分布模型库中的一个模型，对应速分布模型库中所有模型的空域平均波速的最大值，对应速分布模型库中所有模型的空域平均波速的最小值。

8、设定频域拟保留的坐标范围d：即模型库中所有模型的最大d值，此值决定了模型库中每一个模型满足|f(u，v)|≠0的最大频域坐标范围。

9、记模型库中任意一个模型，即一幅空域波速分布为f(x，y)，对其做dft(discretefourier transformation，离散傅里叶变换)得其频域波速分布f(u，v)，记f(u，v)中心点的坐标为(uc，vc)，其中：

10、

11、当d确定时，模型库中所有模型均满足：当或时，|f(u，v)|＝0。

12、设定不同频域坐标处取模上限nmd：此值决定了不同d处|f(u，v)|的上限和|f(u，v)|可取值数量的上限：

13、记坐标(u，v)处的|f(u，v)|为|f(u，v)||u，v，记所有满足max(u，v)＝d的|f(u，v)||u，v为|f(u，v)||d，则|f(u，v)||d可取的值为：

14、

15、其中，(|f(u，v)||u，v)min为所有模型的|f(u，v)||u，v可能取到的最小值，此值满足：

16、(|f(u，v)||u，v)min＝0    (3)

17、(|f(u，v)||u，v)max为所有模型的|f(u，v)||u，v可能取到的最大值，此值满足：

18、

19、式(4)中为最大可能的波速的能量；当某个f(x，y)满足：min(f(x，y))＝max(f(x，y))＝vgmax时，记此f(x，y)为fmax(x，y)，有

20、式(2)中的nm为所有|f(u，v)||u，v可被离散的最大数量且nm为整数；

21、式(2)中的nmd满足：nmd∈[0，nm-1]且nmd为整数。

22、设定模离散数nm和幅角离散数na：

23、nm即式(2)中的nm，是所有|f(u，v)||u，v可离散的最大数量，为整数；

24、na：由dft公式可知，f(u，v)为复数，定义f(u，v)位于坐标(u，v)处的幅角为arg(u，v)，限定arg(u，v)∈[0，2π)，幅角离散数na表示所有arg(u，v)最大可被离散为[0，2π)之间的na个角度，这些角度的最小值为0，相邻两个角度之间的间隔为

25、步骤1-2：根据式(3)～式(4)计算模的最大取值范围：

26、即所有的|f(u，v)||u，v都应满足|f(u，v)||u，v∈[(|f(u，v)||u，v)min，(|f(u，v)||u，v)max]，(|f(u，v)||u，v)min计算公式为式(3)，(|f(u，v)||u，v)max计算公式为式(4)；

27、根据式(2)计算每个模的取值范围：

28、根据nmd和式(2)，计算不同d处的|f(u，v)|所有可能的取值。

29、获得模的所有可能组合的集合m1：

30、由于每个(u，v)处的|f(u，v)|的可能取值均已确定，则所有可能的|f(u，v)|，即集合m1，由各(u，v)处|f(u，v)|的组合获得。

31、步骤1-3：根据式(5)从m1中筛选出集合m2：

32、

33、步骤2：

34、步骤2-1：设定频域坐标(0，0)处能量占比下限emax：

35、定义频域坐标(u，v)处图像的能量比为：

36、

37、emax满足：

38、e(u，v)|u＝0，v＝0≥emax     (7)

39、且有0≤emax≤1。

40、设定能量比阈值et和幅角离散数na：

41、能量比阈值向量et和幅角离散数向量na用于计算不同(u，v)处的幅角的离散数；

42、定义能量比阈值向量

43、et＝[et0，et1，l，eti，etj，l，etnt-1，etnt]    (8)

44、其中et0＝0，etnt＝1，nt≤2d2+2d，且et0＜et1＜l eti＜etj＜l＜etnt；

45、幅角离散数向量

46、na＝[na0，1，na1，2，l，nai，j，l，nant-1，nt]    (9)

47、其中任意元素nai，j表示当eti＜e(u，v)|u，v≤etj时，幅角arg(u，v)|u，v的取值范围被离散为nai，j个值，且有max{na0，1，l，nant-1，nt}≤na，min{na0，1，l，nant-1，nt}≥2。

48、步骤22：根据式(7)从m2中筛选出集合m3：

49、从集合m2中选出所有满足式(7)的模形成集合m3。

50、步骤2-3：根据式(8)～(9)计算m3中每个模对应的幅角的组合：

51、由于et0＝0，etnt＝1，m3中的每个模都能获得与其对应的幅角离散数，将这些幅角离散数相乘，可得所有坐标(u，v)处的幅角离散数，组合这些可能的幅角，可得所有坐标(u，v)处的幅角的可能值。

52、步骤2-4：获得所有模、幅角组合的集合ma1：

53、组合m3中所有可能的模和与每个模对应的所有可能的幅角，即可得集合ma1，ma1为当前参数条件下所有可能的f(u，v)的集合。

54、步骤3：

55、步骤3-1：根据式(10)获得与mal对应的空域图像集合p1。

56、

57、步骤3-2：根据vgmax和vgmin从p1中筛选出空域图像集合p2：

58、即当p1中的f(x，y)满足

59、min(f(x，y))≥vgmin且max(f(x，y))≤vgmax    (11)

60、时，将此f(x，y)选入集合p2。

61、步骤4：根据p2进行正演，获得正演结果集合r：

62、集合p2中保存的是筛选处的空域波速分布f(x，y)，对这些f(x，y)应用二维声波方程可以建立每个f(x，y)对应的震动波场，所有震动波场形成正演结果集合r；

63、上述各集合间的关系：p2为集合p2中元素，满足vgmin≤p2≤vgmax。

64、所述反演过程包括：

65、步骤5：记r中的元素为dtheo，ma2中的元素为mtheo，r中的元素为dtheo，dobs为观测值，反演过程即是根据dobs按照式(12)在r中筛选出dtheo，则与dtheo对应的mtheo构成反演结果集，式(12)如下：

66、

67、为目标函数，ulobj∈(0，1]。

68、进一步的，所述步骤1-1中，目标区域空间的m＝50，n＝50，剖分所用正方形边长为dx＝dy＝5米。

69、进一步的，所述步骤1-2中，获得模的所有可能组合的集合m1：

70、由于每个(u，v)处的|f(u，v)|的可能取值均已确定，则所有可能的|f(u，v)|，即集合m1，由各(u，v)处|f(u，v)|的组合获得；当d＝1时，假设(26，26)处|f(u，v)|有50个，(25，27)、(25，26)、(26，27)、(27，27)处|f(u，v)|有25个，则|f(u，v)|的组合数为50×254＝19531250个，这里当d＝1时仅取了5个坐标点是因为根据dft的公式可知，|f(u，v)|关于远点对称。

71、本发明的有益效果是：本发明基于预计算思想提出一种新型的二维波速分布反演方法。首先基于傅里叶基构建二维波速分布模型库，通过限制目标区域波速均值的最大值与最小值、模型直流分量的能量占比等参数可将模型库大小控制在可接受范围内。然后基于生成的模型库，通过合理设计的目标函数从模型库中筛选出符合条件的所有模型。与传统方法不同的是，本发明所述方法通过检索模型库，给出了其中符合目标函数的所有可能结果，具备更快的反演速度和更高的可靠性，解决了目前主流方法的局限性问题。