CN103675905B - 一种基于优化系数的地震波场模拟方法及装置 - Google Patents

一种基于优化系数的地震波场模拟方法及装置 Download PDF

Info

Publication number
CN103675905B
CN103675905B CN201210343161.1A CN201210343161A CN103675905B CN 103675905 B CN103675905 B CN 103675905B CN 201210343161 A CN201210343161 A CN 201210343161A CN 103675905 B CN103675905 B CN 103675905B
Authority
CN
China
Prior art keywords
interim coefficient
coefficient
current
finite difference
equations
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN201210343161.1A
Other languages
English (en)
Other versions
CN103675905A (zh
Inventor
张金海
姚振兴
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Institute of Geology and Geophysics of CAS
Original Assignee
Institute of Geology and Geophysics of CAS
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Institute of Geology and Geophysics of CAS filed Critical Institute of Geology and Geophysics of CAS
Priority to CN201210343161.1A priority Critical patent/CN103675905B/zh
Priority to US14/117,307 priority patent/US20150134308A1/en
Priority to PCT/CN2012/084083 priority patent/WO2014040338A1/zh
Publication of CN103675905A publication Critical patent/CN103675905A/zh
Application granted granted Critical
Publication of CN103675905B publication Critical patent/CN103675905B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. analysis, for interpretation, for correction
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. analysis, for interpretation, for correction
    • G01V1/30Analysis
    • G01V1/301Analysis for determining seismic cross-sections or geostructures
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. analysis, for interpretation, for correction
    • G01V1/30Analysis
    • G01V1/303Analysis for determining velocity profiles or travel times
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/23Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/60Analysis
    • G01V2210/67Wave propagation modeling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling

Abstract

本发明提供一种优化系数获取方法、装置及相关波场模拟方法、装置,由于本发明通过判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,筛选出符合条件的当前临时系数{Bn}加入待选结果;通过找到当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值均满足第一条件的最大当前离散值确定精度覆盖范围,又将精度覆盖范围最大的一组当前临时系数{Bn}作为第一类优化系数{bn},从而,在随机产生的若干组当前临时系数{Bn}中,查询出一组精度覆盖范围最大的第一类优化系数{bn}作为控制有限差分格式的优化系数,提高了低阶有限差分格式的频率响应范围,使利用优化系数控制的有限差分格式对震源点进行的地震波场模拟效果大大提高。

Description

一种基于优化系数的地震波场模拟方法及装置
技术领域
本发明涉及地球物理勘探领域,特别涉及一种优化系数获取方法、装置及相关波场模拟方法、装置。
背景技术
地震波在各个时间和空间上都是变化的,在野外环境下可通过检波器获得该环境下地震波的实测震动信号,例如,放炮或敲击产生初始激励信号,放炮或敲击位置就是震源点,在地表的一些空间点或者井孔侧壁上放置检波器,获得检波器所在位置的实测震动信号。利用地震波场模拟可以获得野外检波器相同位置的模拟记录,不断的改变地震波传播速度的空间分布,最终使地震波场模拟得到的模拟记录同实测震动信号相一致,实现通过在计算机上模拟震源点周围介质中的波动现象,了解实际地下介质的属性的目的。
可见,地震波场模拟对于与波动现象有关的地震学问题研究具有重要意义,在地震勘探和地震学各工作阶段中都有重要的作用,应用于地震资料采集、处理、解释和地下资源开发工程的各个环节。高精度的地震波场模拟,有助于人们提高复杂勘探目标中地震波传播规律的认识,解决地下矿产资源勘探、开发工作中的各种问题。
地震波场模拟包括以波动方程为基础的地震勘探逆时偏移成像、全波形反演、地震波模拟等等。有限差分法将波动方程中波场函数的空间偏导数和时间偏导数用相应空间和时间的差分来代替,是实现地震波场模拟的主要方法之一,例如:
以波动方程的二阶空间偏导数的有限差分离散为例,对某连续函数f(x)的二阶空间偏导数进行有限差分法离散,实际上是在x=0位置进行如下的泰勒展开:
∂ 2 f ∂ x 2 ≈ 1 Δ 2 Σ n = - N / 2 N / 2 a n N [ - 2 n 2 c o s ( n π ) ] f n
在该式中,偶数N是有限差分格式进行泰勒展开的阶数,Δ是沿空间x方向的空间网格间距,是由如下二项式公式定义的常规系数:
a n N = N N 2 + n / N N 2
由于泰勒展开本身具有局部展开、且收敛速度缓慢的局限性,其主要缺点是对于频带范围较宽的数据有很强的数值频散噪音,而数值频散噪音直接影响了地震波场模拟的精度。实际应用中为了尽可能减少这种噪音的影响,目前主要有两种思路:
(1)采用阶数更高的泰勒展开,即增加更高阶的修正项。参见图1,该图显示的曲线是评价地震波场模拟方法性能的常用手段,横坐标是离散变量波数范围,纵坐标是绝对误差范围;通常认为,绝对误差越小且沿横坐标,即离散变量波数范围跨越的范围越大则代表方法的精度覆盖范围越大,受数值频散误差的影响也越小。如图所示,阶数更高的泰勒展开,其精度越高。但该思路的主要缺点是效果微弱,反而带来成倍的计算量的增加,对数据量大和迭代次数频繁的地震波场模拟,比如地震偏移成像和波形反演等,往往是灾难性的,虽然通过减小空间网格Δ可以有效缓解频散问题,但此时的内存需求又会成倍增加,往往使大规模的三维空间模型很难在现有的计算机条件下进行处理。
(2)直接降低原始信号的主频,即通过滤波消除频率较高的成分来迎合有限差分方法的苛刻要求。这种思路的缺点是直接降低了处理的分辨率,因为高频成分是提高分辨率不可或缺的有效成分。
发明内容
有鉴于此,本发明的主要目的在于提供一种优化系数获取方法、装置及相关波场模拟方法、装置,以实现通过获取精度覆盖范围更广的优化系数代替有限差分的常规系数、利用该优化系数控制有限差分格式提高地震波场模拟的精度的目的。
本发明提供了一种优化系数获取方法,该方法包括:
初始化步骤、计算步骤、检验步骤获取步骤、干扰步骤和输出步骤:
所述初始化步骤包括:
设置误差限T的值;
设置当前离散值的初值;
设置优化系数输出条件;
所述计算步骤包括:
随机产生至少一组当前临时系数{Bn},其中 为Bn预设的浮动上限,为Bn预设的浮动下限,其中所述当前临时系数{Bn}中Bn的个数由有限差分格式具体采用的阶数N决定;
所述检验步骤包括:
判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件;
其中,所述第一条件具体为理想值与实际值之间的差值E小于或者等于预设的误差限T,所述理想值具体为第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C,所述实际值具体为所述第一类方程的空间偏导数在利用当前临时系数Bn控制的有限差分格式的傅里叶变换在离散变量Kx(i)取第i个离散值时的结果,所述离散变量Kx(i)的离散值的范围为0≤Kx(i)<π,C为所述第一类方程的空间偏导数的阶数,为虚数单位;
如果满足第一条件,进入所述获取步骤;
如果不满足第一条件,进入所述干扰步骤;
所述获取步骤包括:
将所述当前临时系数{Bn}加入第一类待选结果;
根据判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取所述当前临时系数{Bn}的精度覆盖范围,所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围的任意离散值均满足第一条件的最大离散值;
所述干扰步骤包括:
判断优化系数输出条件是否满足;
如果优化系数输出条件未满足,将所述当前临时系数{Bn}在当前基础上进行调整,所述当前临时系数{Bn}调整后的值不超过{Bn}预设的浮动上限和下限,更新所述当前临时系数{Bn}为当前临时系数{Bn}调整后的值,进入所述检验步骤;
如果优化系数输出条件满足,进入所述输出步骤;
所述输出步骤包括:
将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}。
优选地,在所述计算步骤中随机产生一组当前临时系数{Bn};
在所述计算步骤之后,进入检验步骤之前,还包括:将当前临时系数{Bn}的值在当前基础上进行调整,调整后的值不超过{Bn}预设的浮动上限和下限,获得调整后临时系数{Bn′};
所述前一临时系数{Bn″}等于当前临时系数{Bn};
所述当前临时系数{Bn}等于所述调整后临时系数{Bn′};
在所述获取步骤中,还包括:所述前一临时系数{Bn″}等于当前临时系数{Bn};
所述初始化步骤还包括:预设温度初值A,预设降温速率α,预设温度最小值A0
所述检验步骤中,如果不满足第一条件,进入所述干扰步骤之前,还包括:判断接受当前解的概率是否大于随机数p,如果否,所述当前临时系数{Bn}等于前一临时系数{Bn″}其中E(当前临时系数)-E(前一临时系数)具体为所述第一类方程的空间偏导数在利用当前临时系数{Bn}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果与所述第一类方程的空间偏导数在利用前一临时系数{Bn″}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果之差,所述随机数p具体为0到1之间的随机数;
所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:判断所述A是否大于A0,如果是,则A=A*α,重新设置优化系数输出条件,重新进入所述干扰步骤;
如果A小于等于A0,进入所述输出步骤。
优选地,所述计算步骤还包括:设置当前离散值为无解状态;
所述获取步骤还包括:设置当前离散值为有解状态;
所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:如果所述A小于等于A0,则判断是否所述当前离散值<π,如果所述当前离散值<π,并且当前离散值为有解状态,则将所述当前离散值增加一个离散间隔作为当前离散值,重新进入所述计算步骤;
如果所述当前离散值≥π或者当前离散值为无解状态,进入所述输出步骤。
优选地,当所述有限差分格式不是交错网格有限差分时,所述预设的误差限T具体为0.0001。
优选地,当所述有限差分是交错网格有限差分时,所述预设的误差限T具体为0.00005。
经过以上本发明提出的获取第一类优化系数{bn}的计算步骤、检验步骤和输出步骤,可以得到以下优选的第一类优化系数:
当所述第一类方程具体为一阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中0.0834≤b-2≤0.1985,-0.1985≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中-0.0357≤b-3≤-0.0167,0.1501≤b-2≤0.2912,-0.2912≤b2≤-0.1501,0.0167≤b3≤0.0357;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中0.0036≤b-4≤0.0097,-0.0669≤b-3≤-0.0381,0.2001≤b-2≤0.3698,-0.3698≤b2≤-0.2001,0.0381≤b3≤0.0669,-0.0097≤b4≤-0.0036;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中-0.0078≤b-5≤-0.0008,0.01≤b-4≤0.0299,-0.1337≤b-3≤-0.0596,0.2381≤b-2≤0.3325,-0.3325≤b2≤-0.2381,0.0596≤b3≤0.1337,-0.0299≤b4≤-0.01,0.0008≤b5≤0.0078;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中0.0001≤b-6≤0.0071,-0.0148≤b-5≤-0.0026,0.0179≤b-4≤0.0588,-0.1527≤b-3≤-0.0794,0.2679≤b-2≤0.3766,-0.3766≤b2≤-0.2679,0.0794≤b3≤0.1527,-0.0588≤b4≤-0.0179,0.0026≤b5≤0.0148,-0.0071≤b6≤-0.0001;
当所述第一类方程具体为一阶偏微分方程,所述有限差分是交错网格有限差分时,
用于控制4阶交错网格的有限差分格式的第一类优化系数bn具体为b-1,b1,b2,其中0.04167≤b-1≤0.0913,-0.0913≤b2≤-0.04167;
用于控制6阶交错网格的有限差分格式的第一类优化系数bn具体为b-2,b-1,b1,b2,b3,其中-0.0761≤b-2≤-0.0047,0.0652≤b-1≤0.1820,-0.1820≤b2≤-0.0652,0.0047≤b3≤0.0761;
用于控制8阶交错网格的有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b1,b2,b3,b4,其中0.0007≤b-3≤0.0034,-0.0188≤b-2≤-0.0096,0.0798≤b-1≤0.1465,-0.1465≤b2≤-0.0798,0.0096≤b3≤0.0188,-0.0034≤b4≤-0.0007;
用于控制10阶交错网格的有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,其中-0.0088≤b-4≤-0.0002,0.0018≤b-3≤0.0084,-0.0139≤b-2≤-0.0298,0.0898≤b-1≤0.1969,-0.1969≤b2≤-0.0898,0.0139≤b3≤0.0298,-0.0084≤b4≤-0.0018,0.0002≤b5≤0.0088;
用于控制12阶交错网格的有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,b6,其中0.0002≤b-5≤0.009,-0.0046≤b-4≤-0.0004,0.0030≤b-3≤0.0979,-0.0599≤b-2≤-0.0175,0.0970≤b-1≤0.1953,-0.1953≤b2≤-0.0970,0.0175≤b3≤0.0599,-0.0979≤b4≤-0.0030,0.0004≤b5≤0.0046,-0.009≤b6≤-0.0002;
当所述第一类方程具体为二阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中-0.1648≤b-2≤-0.0834,-0.1648≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中0.0112≤b-3≤0.0373,-0.3018≤b-2≤-0.1510,-0.3018≤b2≤-0.1510,0.0112≤b3≤0.0373;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中-0.0086≤b-4≤-0.0018,0.0254≤b-3≤0.0585,-0.3855≤b-2≤-0.2001,-0.3855≤b2≤-0.2001,0.0254≤b3≤0.0585,-0.0086≤b4≤-0.0018;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中0.0004≤b-5≤0.0038,-0.0188≤b-4≤-0.0050,0.0397≤b-3≤0.0837,-0.4826≤b-2≤-0.2384,-0.4826≤b2≤-0.2384,0.0397≤b3≤0.0837,-0.0188≤b4≤-0.0050,0.0004≤b5≤0.0038;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中-0.0037≤b-6≤-0.0007,0.0011≤b-5≤0.0077,-0.0327≤b-4≤-0.0090,0.0530≤b-3≤0.1128,-0.3927≤b-2≤-0.2679,-0.3927≤b2≤-0.2679,0.0530≤b3≤0.1128,-0.0327≤b4≤-0.0090,0.0011≤b5≤0.0077,-0.0037≤b6≤-0.0007。
本发明还提供一种优化系数获取装置,该装置包括:
初始化单元:用于设置误差限T的值,设置当前离散值的初值,设置优化系数输出条件;
计算单元:用于随机产生至少一组当前临时系数{Bn},其中 为Bn预设的浮动上限,为Bn预设的浮动下限,其中所述当前临时系数{Bn}中Bn的个数由有限差分格式具体采用的阶数N决定;
检验单元:用于判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件;
其中,所述第一条件具体为理想值与实际值之间的差值E小于或者等于预设的误差限T,所述理想值具体为所述第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C,所述实际值具体为所述第一类方程的空间偏导数在利用当前临时系数Bn控制的有限差分格式的傅里叶变换在离散变量Kx(i)取第i个离散值时的结果,所述离散变量Kx(i)的离散值的范围为0≤Kx(i)<π,C为所述第一类方程的空间偏导数的阶数,为虚数单位;
如果满足第一条件,将所述当前临时系数{Bn}发送至获取单元,触发所述获取单元执行;
如果不满足第一条件,将所述当前临时系数{Bn}发送至干扰单元,触发所述干扰单元执行;
获取单元:用于将所述当前临时系数{Bn}加入第一类待选结果;
根据判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取所述当前临时系数{Bn}的精度覆盖范围,所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;
干扰单元:用于判断优化系数输出条件是否满足;
如果优化系数输出条件未满足,将所述当前临时系数{Bn}在当前基础上进行调整,所述当前临时系数{Bn}调整后的值不超过{Bn}预设的浮动上限和下限,更新所述当前临时系数{Bn}为当前临时系数{Bn}调整后的值,将所述当前临时系数{Bn}发送至所述检验单元,触发所述检验单元执行;
如果优化系数输出条件满足,触发所述输出单元执行;
输出单元:用于将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}。
本发明还提供一种基于优化系数的地震波场模拟方法,该方法包括:
获取震源点激发的波动数据,所述震源点激发的波动数据至少包括震源点波动速度、震源点空间坐标和震源点时间坐标;
获取震源点激发的地震波场模拟涉及的第一类方程;
将所述震源点激发的波动数据作为所述第一类方程的输入数据,应用以上所述一种优化系数获取方法获取的第一类优化系数{bn}控制有限差分格式对震源点激发的地震波场进行模拟。
本发明还提供一种基于优化系数的地震波场模拟装置,该装置包括:
预处理单元:用于获取震源点激发的波动数据,所述震源点激发的波动数据至少包括震源点波动速度、震源点空间坐标和震源点时间坐标;获取震源点激发的地震波场模拟涉及的第一类方程;
模拟单元:用于将所述震源点激发的波动数据作为所述第一类方程的输入数据,应用以上所述一种优化系数获取方法获取的第一类优化系数{bn}控制有限差分格式对震源点激发的地震波场进行模拟。
可见本发明具有如下有益效果:
由于本发明通过判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,筛选出符合条件的当前临时系数{Bn}加入待选结果;通过获取所述当前临时系数{Bn}的精度覆盖范围,找到当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;最后又将精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn},从而,在随机产生的若干组当前临时系数{Bn}中,查询出精度覆盖范围最大的第一类优化系数{bn}作为控制有限差分格式的优化系数,提高了低阶有限差分格式的频率响应范围,使利用优化系数控制的有限差分格式对震源点进行的地震波场模拟效果大大提高;
其次,在当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)的当前离散值不满足第一条件时,又通过模拟退火算法判断接受当前解的概率,跳出该组当前临时系数{Bn}这个局部优化系数,重新计算当前临时系数{Bn},增加搜索到最优化的第一类优化系数{bn}的可能性;
另外,与现有有限差分格式的固定常规系数不同的是,本发明所获取的优化系数{bn}可以通过调整预设误差限T来满足实际应用中不同的精度需求,相对大一些的预设误差限T会使精度覆盖范围得到明显提升,但实际精度将会比阀值小的略低,因此可以根据具体应用的实际需要来合理选择;
而且经过实验数据验证,在效果相当的前提下,依据本发明获取的第一类优化系数{bn},其控制有限差分格式进行地震波场模拟,在耗费内存和计算量方面,较常规的有限差分方法有明显的降低。
附图说明
图1是现有常规系数控制有限差分格式的精度覆盖范围图例;
图2是本发明一种优化系数获取方法步骤图例;
图3是本发明一种优化系数获取方法优选实施例的步骤图例;
图4是本发明一种优化系数获取装置的组成图例;
图5是本发明一种基于优化系数的地震波场模拟方法步骤图例;
图6是本发明一种基于优化系数的地震波场模拟装置组成图例;
图7-1是本发明实验一,现有常规系数控制有限差分格式的精度覆盖范围图例;
图7-2是本发明实验一,优化系数控制有限差分格式的精度覆盖范围图例;
图8-1是本发明实验二,地震波场模拟采取Marmousi模型模拟效果图例;
图8-2是本发明实验二,对地震波场模拟采取Marmousi模型,采取现有常规系数控制有限差分格式与本发明优化系数控制有限差分格式的精度随时间变化曲线图例;
图9-1是本发明实验三,现有常规系数控制有限差分格式进行地震波场模拟的内存耗费量,计算量比例示意图例;
图9-2是本发明实验三,本发明优化系数控制有限差分格式进行地震波场模拟的内存耗费量,计算量比例示意图例。
具体实施方式
为使本发明的上述目的、特征和优点能够更加明显易懂,下面结合附图和具体实施方式对本发明实施例作进一步详细的说明。
本发明提供了一种优化系数获取方法,参见图2,该方法包括:初始化步骤、计算步骤、检验步骤、获取步骤、干扰步骤和输出步骤:
S201、所述初始化步骤包括:
设置误差限T的值;
设置当前离散值的初值;
设置优化系数输出条件;
S202、所述计算步骤包括:
S202.1、随机产生至少一组当前临时系数{Bn},其中 为Bn预设的浮动上限,为Bn预设的浮动下限,其中所述当前临时系数{Bn}中Bn的个数由有限差分格式具体采用的阶数N决定;
S203、所述检验步骤包括:
S203.1、判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件;
其中,所述第一条件具体为理想值与实际值之间的差值E小于或者等于预设的误差限T,所述理想值具体为所述第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C,所述实际值具体为所述第一类方程的空间偏导数在利用当前临时系数Bn控制的有限差分格式的傅里叶变换在离散变量Kx(i)取第i个离散值时的结果,所述离散变量Kx(i)的离散值的范围为0≤Kx(i)<π,C为所述第一类方程的空间偏导数的阶数,为虚数单位;
如果满足第一条件,进入所述获取步骤;
如果不满足第一条件,进入所述干扰步骤;
需要说明的是,本发明所述离散变量以波动方程为例,离散变量具体为离散的波数,范围为π,离散变量的离散值之间的离散间隔,应该预设为一个相对较小的间隔,如该离散间隔越小,计算量也随之增大,在本发明的一个实施例中,所述离散间隔具体为
S204、所述获取步骤包括:
S204.1、将所述当前临时系数{Bn}加入第一类待选结果;
S204.2、根据判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取所述当前临时系数{Bn}的精度覆盖范围,所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;
S205、所述干扰步骤包括:
S205.1、判断优化系数输出条件是否满足,如果优化系数输出条件满足,进入所述输出步骤;
S205.2、如果优化系数输出条件未满足,将所述当前临时系数{Bn}在当前基础上进行调整,所述当前临时系数{Bn}调整后的值不超过{Bn}预设的浮动上限和下限,更新所述当前临时系数{Bn}为当前临时系数{Bn}调整后的值,进入所述检验步骤;
S206、所述输出步骤包括:
将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}。
需要说明的是,以上各个步骤中所述将当前临时系数{Bn}在当前基础上进行调整,具体可以按照需要或者经验进行调整,例如,以下三种实施方式:
(1)按照所述计算步骤的方法对当前临时系数{Bn}进行随机运算;
(2)预设固定的浮动百分比,将当前临时系数{Bn}在当前基础上浮动一定百分比,例如上下浮动10%;
(3)预设随当前临时系数{Bn}调整次数变化的浮动百分比,例如,第一次调整,上下浮动20%,第二次,上下浮动19.5%,第三次,上下浮动19%,第四次,上下浮动18.5%,依次类推逐次降低浮动大小,以达到搜索范围逐渐收敛的效果。
通过以上S201到S206步骤可见,由于本发明通过判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,筛选出符合条件的当前临时系数{Bn}加入待选结果;通过获取所述当前临时系数{Bn}的精度覆盖范围,找到当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;最后又将精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn},从而,在随机产生的若干组当前临时系数{Bn}中,查询出精度覆盖范围最大的第一类优化系数{bn}作为控制有限差分格式的优化系数,提高了低阶有限差分格式的频率响应范围,使利用优化系数控制的有限差分格式对震源点进行的地震波场模拟效果大大提高。
下面对步骤S204.2中所述根据判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取当前临时系数{Bn}的精度覆盖范围进行说明:
所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;
根据本发明所述根据判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件的提示,获取当前临时系数{Bn}的精度覆盖范围,可以通过以下赋值步骤和筛选步骤获取所述当前临时系数{Bn}的精度覆盖范围:
所述赋值步骤包括:
设置当前临时离散值等于所述当前离散值;
所述筛选步骤包括:
判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前临时离散值,是否均满足第一条件;
若判断出当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前临时离散值,均满足第一条件,将当前临时离散值作为前一离散值,将当前临时离散值增加一个离散间隔作为当前离散值,重新进入所述筛选步骤;
若判断出当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)在当前临时离散值不满足第一条件,将所述前一离散值作为所述当前临时系数{Bn}的精度覆盖范围。
需要说明的是,为了尽量增加搜索到最优化的第一类优化系数{bn}的可能性,本发明还提出利用模拟退火算法进一步搜索第一类优化系数{bn},为了更清楚的说明该优选实施例的实现过程,下面将该优选实施例获取所述第一类优化系数{bn}的初始化步骤、计算步骤、检验步骤、获取步骤、干扰步骤和输出步骤进行整体的详细说明,参见图3:
S301、所述初始化步骤包括:
设置误差限T的值;
设置当前离散值的初值;
设置优化系数输出条件;
设置温度初值A;
设置降温速率α;
设置温度最小值A0
S302、所述计算步骤包括:
S302.1、随机产生一组当前临时系数{Bn},其中 为Bn预设的浮动上限,为Bn预设的浮动下限,其中所述当前临时系数{Bn}中Bn的个数由有限差分格式具体采用的阶数N决定;
S302.2、将当前临时系数{Bn}的值在当前基础上进行调整,调整后的值不超过{Bn}预设的浮动上限和下限,获得调整后临时系数{Bn′};
所述前一临时系数{Bn″}等于当前临时系数{Bn};
所述当前临时系数{Bn}等于所述调整后临时系数{Bn′};
S303、所述检验步骤包括:
S303.1、判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,如果满足第一条件,进入所述获取步骤;
其中所述第一条件与以上其他实施例中所述含义相同,在此不再赘述;
S303.2、如果不满足第一条件,判断接受当前解的概率是否大于随机数p;
其中E(当前临时系数)-E(前一临时系数)具体为所述第一类方程的空间偏导数在利用当前临时系数{Bn}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果与所述第一类方程的空间偏导数在利用前一临时系数{Bn″}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果之差,所述随机数p具体为0到1之间的随机数;
S303.3、如果否,所述当前临时系数{Bn}等于前一临时系数{Bn″};
S303.4、进入所述干扰步骤;
S304、所述获取步骤包括:
S304.1、将所述当前临时系数{Bn}加入第一类待选结果;
S304.2、根据判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取所述当前临时系数{Bn}的精度覆盖范围,所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;
S304.3、所述前一临时系数{Bn″}等于当前临时系数{Bn};
S305、所述干扰步骤包括:
S305.1、判断优化系数输出条件是否满足;
S305.2、如果优化系数输出条件满足,判断所述A是否大于A0,如果A小于等于A0,进入所述输出步骤;
S305.2a、如果A大于A0,则A=A*α,重新设置优化系数输出条件,重新进入所述干扰步骤;
S305.3、如果优化系数输出条件未满足,将所述当前临时系数{Bn}在当前基础上进行调整,所述当前临时系数{Bn}调整后的值不超过{Bn}预设的浮动上限和下限,更新所述当前临时系数{Bn}为当前临时系数{Bn}调整后的值,进入所述检验步骤;
S306、所述输出步骤包括:
S306.1、将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}。
通过对以上获取所述第一类优化系数{bn}的优选实施例的初始化步骤、计算步骤、检验步骤、获取步骤、干扰步骤和输出步骤的详细说明,可见该优选实施例相对于本文上一获取第一类优化系数{bn}的实施例的不同之处在于:
(1)在所述计算步骤中随机产生一组当前临时系数{Bn};
(2)在所述计算步骤之后,进入检验步骤之前,还包括:将当前临时系数{Bn}的值在当前基础上进行调整,调整后的值不超过{Bn}预设的浮动上限和下限,获得调整后临时系数{Bn′};
所述前一临时系数{Bn″}等于当前临时系数{Bn};
所述当前临时系数{Bn}等于所述调整后临时系数{Bn′};
(3)在所述获取步骤中,还包括:所述前一临时系数{Bn″}等于当前临时系数{Bn};
(4)所述初始化步骤还包括:预设温度初值A,预设降温速率α,预设温度最小值A0
(5)所述检验步骤中,如果不满足第一条件,进入所述干扰步骤之前,还包括:判断接受当前解的概率是否大于随机数p,如果否,所述当前临时系数{Bn}等于前一临时系数{Bn″},其中E(当前临时系数)-E(前一临时系数)具体为所述第一类方程的空间偏导数在利用当前临时系数{Bn}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果与所述第一类方程的空间偏导数在利用前一临时系数{Bn″}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果之差,所述随机数p具体为0到1之间的随机数;
(6)所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:判断所述A是否大于A0,如果是,则A=A*α,重新设置优化系数输出条件,重新进入所述干扰步骤;
如果A小于等于A0,进入所述输出步骤。
而且,在本发明又一实施例中,提出在上述优选实施例的基础之上,通过以下步骤实现进一步扩大搜索到最优化的第一类优化系数{bn}的可能,包括:
所述计算步骤还包括:设置当前离散值为无解状态;
所述获取步骤还包括:设置当前离散值为有解状态,判断是否所述当前离散值<π,如果是,将所述当前离散值增加一个离散间隔作为当前离散值,重新进入所述计算步骤,如果否,进入所述输出步骤;
所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:如果所述A小于等于A0,则判断是否所述当前离散值<π,如果所述当前离散值<π,并且当前离散值为有解状态,则将所述当前离散值增加一个离散间隔作为当前离散值,重新进入所述计算步骤;
如果所述当前离散值≥π或者当前离散值为无解状态,进入所述输出步骤。
可见,该优选实施例相当于在模拟退火算法降温流程外层又增加了一层寻找第一类优化系数{bn}的循环,通过判断当前离散值内是否有解,决定是否在下一离散值内继续寻找第一类优化系数{bn};在当前离散值之内有解的情况下,通过逐步增加当前离散值,扩大寻找到最优化的第一类优化系数{bn}可能。
需要说明的是,对于所述输出步骤将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn},若第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}为多组,为了选择出最优化的一组当前临时系数{Bn},本发明还包括:
通过计算所述第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C与所述第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的结果之差,获得第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的误差。
在计算出当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的误差基础之上,所述输出步骤将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}具体为:
将第一类待选结果中精度覆盖范围最大,且误差和最小的当前临时系数{Bn}作为第一类优化系数{bn};
所述当前临时系数{Bn}的误差和通过计算所述第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的误差之和获得。
下面,对本发明该优选实施例中所述优化系数输出条件通过以下两个实施例进行详细说明:
在本发明的一个实施例中,该步骤中所述优化系数输出条件具体可以为重新进入计算步骤的次数超过预设的干扰次数阀值,可见,该干扰次数阀值越大,跳出当前临时系数{Bn}局部优化系数,重新计算当前临时系数{Bn}的次数就越多,增加搜索到最优化的第一类优化系数{bn}的可能性就越大,因此该干扰次数阀值应该预设为一个尽量大的数值,例如,干扰次数阀值预设为60000。
在本发明的又一实施例中,该步骤中所述优化系数输出条件具体可以为所述第一类优化系数{bn}的获取时间超过预设的时间阀值,同样可见,该时间阀值越大,跳出当前临时系数{Bn}局部优化系数,重新计算当前临时系数{Bn}的机会就越多,增加搜索到最优化的第一类优化系数{bn}的可能性就越大,因此该时间阀值应该预设为一个尽量大的数值,例如,时间阀值预设为7天,通过该实施例,可以使本发明的地震波场模拟实现的效果与实际工作时间上的需求相适应。
考虑到在本发明具体实施到不同应用场景中,对最终获取的第一类优化系数{bn}精度覆盖范围的最低要求是不同的,因此,在本发明的一个优选实施例中,还包括:
预设精度覆盖范围阀值;
在进入所述输出步骤之前,还包括:
判断在第一类待选结果中,精度覆盖范围最大的一组当前临时系数{Bn}的精度覆盖范围是否小于预设的精度覆盖范围阀值,如果是,则放宽Bn预设的浮动上限和/或放宽Bn预设的浮动下限重新进入计算步骤。
当然,放宽Bn预设的浮动上限和/或放宽Bn预设的浮动下限应该参考现有控制有限差分格式的常规系数上下浮动一定范围来进行预设以提高第一类优化系数{bn}获取的效率,一般放宽Bn预设的浮动上限和/或放宽Bn预设的浮动下限为原控制有限差分格式的常规系数的20%~30%。
对于步骤S202.1,在本发明的一个优选实施例中,所述随机产生当前临时系数{Bn}还通过限定所述第一类优化系数{bn}满足一定的优化条件来提高的计算速度和精度,下面,按照所述第一类方程和有限差分格式的类型分为三种情况,详细描述该优选实施例:
(一)在所述第一类方程为一阶偏微分方程,所述有限差分格式不是交错网格有限差分时,本发明还包括限定第一类优化系数{bn}需要满足的优化条件,具体为:
(1)限定所述当前临时系数{Bn}包括第一类临时系数{B-m}、中间临时系数B0和第二类临时系数{Bm},其中m>0;
例如:有限差分格式具体采用6阶,所述当前临时系数{Bn}包括B-3,B-2,B-1,B0,B1,B2,B3
(2)限定所述第一类临时系数{B-m}和第二类临时系数{Bm}以中间临时系数B0为中心奇对称;
(3)限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
根据以上优化条件,以6阶有限差分格式为例,B-3=-B3,B-2=-B2,B-1=-B1
(4)限定所述当前临时系数{Bn}的总和为0;
根据该优化条件,B0=0;
(5)限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,越邻近中间临时系数B0的系数的绝对值越大。
(二)在所述第一类方程为二阶偏微分方程,所述有限差分格式不是交错网格有限差分时,本发明还包括限定第一类优化系数{bn}需要满足的优化条件,具体为:
(1)限定所述当前临时系数{Bn}包括第一类临时系数{B-m}、中间临时系数B0和第二类临时系数{Bm},其中m>0;
例如:有限差分格式具体采用6阶,所述当前临时系数{Bn}包括B-3,B-2,B-1,B0,B1,B2,B3
(2)限定所述第一类临时系数{B-m}和第二类临时系数{Bm}以中间临时系数B0为中心偶对称;
依据该优化条件,以6阶有限差分格式为例,B-3=B3,B-2=B2,B-1=B1
(3)限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
(4)限定所述当前临时系数{Bn}的总和为0;
(5)限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,越邻近中间临时系数B0的系数的绝对值越大。
对于以上(一)、(二)两种情况,所述计算步骤的随机产生至少一组满足所述优化条件的当前临时系数{Bn},具体通过以下步骤实现:
对应每个待求的第二类临时系数Bm各分配一个第一类随机数rm,其中0≤rm≤1;
根据计算出所述第二类临时系数{Bm}的值,其中为Bm预设的浮动上限,为Bm预设的浮动下限;
根据第一类优化系数{bn}的优化条件和所述第二类临时系数{Bm}的值,求出第一类临时系数{B-m}与中间临时系数B0的值。
(三)在所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分时,本发明还包括限定第一类优化系数{bn}需要满足的优化条件,具体为:
(1)限定所述当前临时系数{Bn}包括第一类临时系数{B-m+1}、中间临时系数B1和第二类临时系数{Bm},其中m>1;
例如:有限差分格式具体采用6阶,所述当前临时系数{Bn}包括B-2,B-1,B1,B2,B3
(2)限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}以中间临时系数B1为中心奇对称;
依据该优化条件,以6阶有限差分格式为例,B-2=B3,B-1=B2
(3)限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
(4)限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}中,越邻近中间临时系数B1的系数的绝对值越大。
对于以上(三)的情况,所述计算步骤的随机产生至少一组满足所述优化条件的当前临时系数{Bn},具体通过以下步骤实现:
对应每个待求的第二类临时系数Bm各分配一个第一类随机数rm,其中0≤rm≤1;
根据计算出所述第二类临时系数{Bm}的值,其中为Bm预设的浮动上限,为Bm预设的浮动下限;
根据第一类优化系数{bn}的优化条件和所述第二类临时系数{Bm}的值,求出第一类临时系数{B-m+1}与中间临时系数B1的值。
需要说明的是,在所述计算步骤中,预设的浮动上限和预设的浮动下限可以参考现有控制有限差分格式的常规系数上下浮动一定范围来进行预设。
下面,对本发明步骤S201中,所述预设的误差限T进行详细说明:
在本发明的一个实施例中,有限差分格式不是交错网格有限差分,所述预设的误差限T具体可以为0.0001,在本发明的又一实施例中,有限差分格式是交错网格有限差分,所述预设的误差限T具体可以为0.00005,当然,以上建议的预设的误差限取值附近的某个小数也是可以选择的对象,具体依据本发明具体实施的需求进行设置。但是,合理的选择误差限至关重要,对于过于小的误差限,将会导致精度覆盖的波数范围有限,对于过大的误差限,虽然可以很容易的使精度覆盖的波数范围较大,但是,会对实际应用带来潜在危害,例如,我们的实验结果表明:选取0.0003~0.03的误差限,最大可以覆盖全部的波数范围,但以此误差范围为约束得到的优化系数,其实际的精度较低。因此,不能单纯的通过放大误差限来扩大波数覆盖范围。经过数值实验和理论分析,保证精度,又保证精度覆盖的波数范围较大的误差限为以上所建议的值,即:有限差分格式不是交错网格有限差分,建议预设的误差限T为0.0001,有限差分格式是交错网格有限差分,建议预设的误差限T为0.00005。
下面,再对对本发明步骤S203.1中,所述判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,这一步骤通过本发明的以下三个实施例进行说明:
(一)在本发明所述第一类方程为一阶偏微分方程,所述有限差分格式不是交错网格有限差分这个实施例中:
对某连续函数f(x)的一阶空间偏导数进行有限差分法离散,实际上是在x=0位置进行如下格式的泰勒展开:
∂ f ∂ x ≈ 1 Δ Σ n = - N / 2 N / 2 b n c o s ( n π ) f n
因此,所述判断利用当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
其中Δ为震源点激发的波动数据的空间网格间距;
(二)在本发明所述第一类方程为二阶偏微分方程,所述有限差分格式不是交错网格有限差分这个实施例中,
对某连续函数f(x)的二阶空间偏导数进行有限差分法离散,实际上是在x=0位置进行如下格式的泰勒展开:
∂ 2 f ∂ x 2 ≈ 1 Δ 2 Σ n = - N / 2 N / 2 b n c o s ( n π ) f n
因此,所述判断利用当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( K x ( i ) , T ) ≡ m a x 0 ≤ k x ( i ) | - K x ( i ) 2 Δ 2 - Σ n = - N / 2 N / 2 B n c o s ( nK x ( i ) Δ ) | ≤ T ;
(三)在本发明所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分这个实施例中,
对某连续函数f(x)的一阶空间偏导数进行有限差分法离散,实际上是在x=0位置进行如下格式的泰勒展开:
∂ f ∂ x ≈ 1 Δ Σ n = - N / 2 N / 2 b n s i n [ ( 0.5 - n ) π ] f n
因此,所述判断利用当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( k x ( i ) , T ) ≡ max 0 ≤ k x ( i ) | - k x ( i ) Δ - Σ n = - N / 2 N / 2 b n sin [ ( 0.5 - n ) k x ( i ) Δ ] | ≤ T .
本发明通过实验数据验证,经过以上获取第一类优化系数{bn}的初始化步骤、计算步骤、检验步骤、获取步骤、干扰步骤和输出步骤,可以得到以下范围的第一类优化系数{bn},使得地震波场模拟效果大大提升:
(一)所述第一类方程具体为一阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中0.0834≤b-2≤0.1985,-0.1985≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中-0.0357≤b-3≤-0.0167,0.1501≤b-2≤0.2912,-0.2912≤b2≤-0.1501,0.0167≤b3≤0.0357;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中0.0036≤b-4≤0.0097,-0.0669≤b-3≤-0.0381,0.2001≤b-2≤0.3698,-0.3698≤b2≤-0.2001,0.0381≤b3≤0.0669,-0.0097≤b4≤-0.0036;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中-0.0078≤b-5≤-0.0008,0.01≤b-4≤0.0299,-0.1337≤b-3≤-0.0596,0.2381≤b-2≤0.3325,-0.3325≤b2≤-0.2381,0.0596≤b3≤0.1337,-0.0299≤b4≤-0.01,0.0008≤b5≤0.0078;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中0.0001≤b-6≤0.0071,-0.0148≤b-5≤-0.0026,0.0179≤b-4≤0.0588,-0.1527≤b-3≤-0.0794,0.2679≤b-2≤0.3766,-0.3766≤b2≤-0.2679,0.0794≤b3≤0.1527,-0.0588≤b4≤-0.0179,0.0026≤b5≤0.0148,-0.0071≤b6≤-0.0001;
(二)所述第一类方程具体为一阶偏微分方程,所述有限差分是交错网格有限差分时,
用于控制4阶交错网格的有限差分格式的第一类优化系数bn具体为b-1,b1,b2,其中0.04167≤b-1≤0.0913,-0.0913≤b2≤-0.04167;
用于控制6阶交错网格的有限差分格式的第一类优化系数bn具体为b-2,b-1,b1,b2,b3,其中-0.0761≤b-2≤-0.0047,0.0652≤b-1≤0.1820,-0.1820≤b2≤-0.0652,0.0047≤b3≤0.0761;
用于控制8阶交错网格的有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b1,b2,b3,b4,其中0.0007≤b-3≤0.0034,-0.0188≤b-2≤-0.0096,0.0798≤b-1≤0.1465,-0.1465≤b2≤-0.0798,0.0096≤b3≤0.0188,-0.0034≤b4≤-0.0007;
用于控制10阶交错网格的有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,其中-0.0088≤b-4≤-0.0002,0.0018≤b-3≤0.0084,-0.0139≤b-2≤-0.0298,0.0898≤b-1≤0.1969,-0.1969≤b2≤-0.0898,0.0139≤b3≤0.0298,-0.0084≤b4≤-0.0018,0.0002≤b5≤0.0088;
用于控制12阶交错网格的有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,b6,其中0.0002≤b-5≤0.009,-0.0046≤b-4≤-0.0004,0.0030≤b-3≤0.0979,-0.0599≤b-2≤-0.0175,0.0970≤b-1≤0.1953,-0.1953≤b2≤-0.0970,0.0175≤b3≤0.0599,-0.0979≤b4≤-0.0030,0.0004≤b5≤0.0046,-0.009≤b6≤-0.0002;
(三)所述第一类方程具体为二阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中-0.1648≤b-2≤-0.0834,-0.1648≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中0.0112≤b-3≤0.0373,-0.3018≤b-2≤-0.1510,-0.3018≤b2≤-0.1510,0.0112≤b3≤0.0373;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中-0.0086≤b-4≤-0.0018,0.0254≤b-3≤0.0585,-0.3855≤b-2≤-0.2001,-0.3855≤b2≤-0.2001,0.0254≤b3≤0.0585,-0.0086≤b4≤-0.0018;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中0.0004≤b-5≤0.0038,-0.0188≤b-4≤-0.0050,0.0397≤b-3≤0.0837,-0.4826≤b-2≤-0.2384,-0.4826≤b2≤-0.2384,0.0397≤b3≤0.0837,-0.0188≤b4≤-0.0050,0.0004≤b5≤0.0038;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中-0.0037≤b-6≤-0.0007,0.0011≤b-5≤0.0077,-0.0327≤b-4≤-0.0090,0.0530≤b-3≤0.1128,-0.3927≤b-2≤-0.2679,-0.3927≤b2≤-0.2679,0.0530≤b3≤0.1128,-0.0327≤b4≤-0.0090,0.0011≤b5≤0.0077,-0.0037≤b6≤-0.0007。
本发明还提供一种优化系数获取装置,参见图4,该装置包括:
初始化单元401:用于设置误差限T的值,设置当前离散值的初值,设置优化系数输出条件;
计算单元402:用于随机产生至少一组当前临时系数{Bn},其中 为Bn预设的浮动上限,为Bn预设的浮动下限,其中所述当前临时系数{Bn}中Bn的个数由有限差分格式具体采用的阶数N决定;
检验单元403:用于判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件;
其中,所述第一条件具体为理想值与实际值之间的差值E小于或者等于预设的误差限T,所述理想值具体为所述第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C,所述实际值具体为所述第一类方程的空间偏导数在利用当前临时系数Bn控制的有限差分格式的傅里叶变换在离散变量Kx(i)取第i个离散值时的结果,所述离散变量Kx(i)的离散值的范围为0≤Kx(i)<π,C为所述第一类方程的空间偏导数的阶数,为虚数单位;
如果满足第一条件,将所述当前临时系数{Bn}发送至获取单元,触发所述获取单元404执行;
如果不满足第一条件,将所述当前临时系数{Bn}发送至干扰单元,触发所述干扰单元405执行;
获取单元404:用于将所述当前临时系数{Bn}加入第一类待选结果;
根据判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取所述当前临时系数{Bn}的精度覆盖范围,所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;
干扰单元405:用于判断优化系数输出条件是否满足;
如果优化系数输出条件未满足,将所述当前临时系数{Bn}在当前基础上进行调整,所述当前临时系数{Bn}调整后的值不超过{Bn}预设的浮动上限和下限,更新所述当前临时系数{Bn}为当前临时系数{Bn}调整后的值,将所述当前临时系数{Bn}发送至所述检验单元403,触发所述检验单元403执行;
如果优化系数输出条件满足,触发所述输出单元406执行;
输出单元406:用于将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}。
本发明还提供一种基于优化系数的地震波场模拟方法,参见图5,该方法包括:
S501、获取震源点激发的波动数据,所述震源点激发的波动数据至少包括震源点波动速度、震源点空间坐标和震源点时间坐标;
S502、获取震源点激发的地震波场模拟涉及的第一类方程;
S503、将所述震源点激发的波动数据作为所述第一类方程的输入数据,应用如以上各实施例所述一种优化系数获取方法各实施例获取的第一类优化系数{bn}控制有限差分格式对震源点激发的地震波场进行模拟。
本发明还提供一种基于优化系数的地震波场模拟装置,参见图6,该装置包括:
预处理单元601:用于获取震源点激发的波动数据,所述震源点激发的波动数据至少包括震源点波动速度、震源点空间坐标和震源点时间坐标;获取震源点激发的地震波场模拟涉及的第一类方程;
模拟单元602:用于将所述震源点激发的波动数据作为所述第一类方程的输入数据,应用如以上所述一种优化系数获取方法各实施例获取的第一类优化系数{bn}控制有限差分格式对震源点激发的地震波场进行模拟。
为了进一步说明本发明的有益效果,通过以下实验数据图例对本发明获取的优化系数{bn}控制的有限差分格式与现有常规系数控制的有限差分格式进行地震波场模拟效果的比较:
(实验一)从精度覆盖波数范围来看:
见图7-1,在图7-1中,横坐标为离散变量波数的范围,纵坐标为有限差分格式的精度,坐标系中的实线曲线为现有常规系数控制的有限差分格式分别在4阶、8阶、12阶、16阶、20阶、24阶、28阶泰勒展开时精度与其覆盖波数的曲线。
见图7-2,在图7-2中,横坐标为离散变量波数的范围,纵坐标为有限差分格式的精度,坐标系中的虚线曲线为本发明获取的优化系数{bn}控制的有限差分格式分别在4阶,8阶泰勒展开时精度与其覆盖波数的曲线。
比较图7-1和图7-2可见,采用本发明获取的优化系数{bn}控制有限差分格式,同现有常规系数控制的有限差分格式相比,在相同阶数的泰勒展开下,本发明获取的优化系数{bn}控制的有限差分格式具有更大的精度覆盖范围,比如,优化系数{bn}控制的有限差分格式的8阶泰勒展开的精度覆盖范围同现有常规系数控制的有限差分格式的12阶泰勒展开的精度覆盖范围基本一致;优化系数{bn}控制的有限差分格式的12阶泰勒展开的精度覆盖范围同现有常规系数控制的有限差分格式的24阶泰勒展开的精度覆盖范围基本一致。
(实验二)从地震波场模拟的精度随着时间变化来看:
见图8-1,该地震波场模拟采取Marmousi模型进行地震波场模拟,为了便于比较,将Marmousi模型的网格统一设定为均匀网格,空间网格间距Δ=5米,模型网格为737×751,Ricker子波的主频为50赫兹,震源点位于水平2000米,纵深20米处,接收点位于水平3000米,纵深5米处;
见图8-2,在图8-2中:
横坐标为时间范围;
纵坐标为地震波场模拟精度范围;
虚线曲线为采用常规系数控制有限差分格式在36阶泰勒展开时,进行地震波场模拟的精度变化曲线;
实线1为现有常规系数控制有限差分格式在12阶泰勒展开时,进行地震波场模拟的精度变化曲线;
实线2为现有常规系数控制有限差分格式在24阶泰勒展开时,进行地震波场模拟的精度变化曲线;
实线3为本发明获取的优化系数{bn}控制的有限差分格式在12阶泰勒展开时,进行地震波场模拟的精度变化曲线;
需要说明的是,图8-2是评价地震波场模拟方法性能的常用手段,将现有常规系数的有限差分格式在36阶泰勒展开进行地震波场模拟的精度曲线,即图8-2中虚线,作为理想值参考,若实线越与虚线一致,说明精度越高;
从图8-2可见,本发明获取的优化系数{bn}控制的有限差分格式在12阶泰勒展开时,进行地震波场模拟的精度变化曲线远远优于现有常规系数控制的有限差分格式在12阶泰勒展开时进行地震波场模拟的精度变化曲线,而且几乎和现有常规系数控制的有限差分格式在24阶泰勒展开时进行地震波场模拟的精度变化曲线相当。
(实验三)从耗费内存和计算量来看:
在该实验中,作为比较的前提,对于给定的速度模型进行地震波场模拟,其尺寸是固定的,但划分网格的间距和网格数目是可以变化的,模型在划分过程中确保不出现数值频散;
在柱状图9-1中,横坐标为有限差分格式泰勒展开的阶数,纵坐标为内存量或者计算量的比例,实心柱1为现有常规系数控制的有限差分格式进行地震波场模拟的内存耗费量,空心柱1为现有常规系数控制的有限差分格式进行地震波场模拟的计算量;
在柱状图9-2中,横坐标为有限差分格式泰勒展开的阶数,,纵坐标为内存量或者计算量的比例,实线柱2为本发明获取的优化系数{bn}控制的有限差分格式进行地震波场模拟的内存耗费量,虚线柱2为本发明获取的优化系数{bn}控制的有限差分格式进行地震波场模拟的计算量;
比较图9-1和图9-2可见,
本发明获取的优化系数{bn}控制的有限差分格式在8阶泰勒展开时地震波场模拟耗费的内存量和计算量,与现有常规系数控制的有限差分格式在12阶泰勒展开时地震波场模拟耗费的内存量和计算量相比较,地震波场模拟耗费的内存量相当,计算量小;
本发明获取的优化系数{bn}控制的有限差分格式在12阶泰勒展开时地震波场模拟耗费的内存量和计算量,与现有常规系数控制的有限差分格式在24阶泰勒展开时地震波场模拟耗费的内存量和计算量相比较,地震波场模拟耗费的内存量相当,计算量小。
需要说明的是,在本文中,诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来,而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且,术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含,从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下,由语句“包括一个……”限定的要素,并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。
以上所述仅为本发明的较佳实施例而已,并非用于限定本发明的保护范围。凡在本发明的精神和原则之内所作的任何修改、等同替换、改进等,均包含在本发明的保护范围内。

Claims (24)

1.一种基于优化系数的地震波场模拟方法,其特征在于,包括:
获取震源点激发的波动数据,所述震源点激发的波动数据至少包括模型介质的波动速度、震源点空间坐标和震源点时间坐标;
获取震源点激发的地震波场模拟涉及的第一类方程;
将所述震源点激发的波动数据作为所述第一类方程的输入数据,应用一种优化系数获取方法获取的第一类优化系数{bn}控制有限差分格式对震源点激发的地震波场进行模拟;
所述优化系数获取方法包括初始化步骤、计算步骤、检验步骤、获取步骤、干扰步骤和输出步骤:
所述初始化步骤包括:
设置误差限T的值;
设置当前离散值的初值;
设置优化系数输出条件;
所述计算步骤包括:
随机产生至少一组当前临时系数{Bn},其中 为Bn预设的浮动上限,为Bn预设的浮动下限,其中所述当前临时系数{Bn}中Bn的个数由有限差分格式具体采用的阶数N决定;
所述检验步骤包括:
判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件;
其中,所述第一条件具体为理想值与实际值之间的差值E小于或者等于预设的误差限T,所述理想值具体为第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C,所述实际值具体为所述第一类方程的空间偏导数在利用当前临时系数Bn控制的有限差分格式的傅里叶变换在离散变量Kx(i)取第i个离散值时的结果,所述离散变量Kx(i)的离散值的范围为0≤Kx(i)<π,C为所述第一类方程的空间偏导数的阶数,为虚数单位;
如果满足第一条件,进入所述获取步骤;
如果不满足第一条件,进入所述干扰步骤;
所述获取步骤包括:
将所述当前临时系数{Bn}加入第一类待选结果;
根据判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取所述当前临时系数{Bn}的精度覆盖范围,所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;
所述干扰步骤包括:
判断优化系数输出条件是否满足;
如果优化系数输出条件未满足,将所述当前临时系数{Bn}在当前基础上进行调整,所述当前临时系数{Bn}调整后的值不超过{Bn}预设的浮动上限和下限,更新所述当前临时系数{Bn}为当前临时系数{Bn}调整后的值,进入所述检验步骤;
如果优化系数输出条件满足,进入所述输出步骤;
所述输出步骤包括:
将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}。
2.根据权利要求1所述的方法,其特征在于,
在所述计算步骤中随机产生一组当前临时系数{Bn};
在所述计算步骤之后,进入检验步骤之前,还包括:将当前临时系数{Bn}的值在当前基础上进行调整,调整后的值不超过{Bn}预设的浮动上限和下限,获得调整后临时系数{Bn′};
前一临时系数{Bn″}等于当前临时系数{Bn};
所述当前临时系数{Bn}等于所述调整后临时系数{Bn′};
在所述获取步骤中,还包括:所述前一临时系数{Bn″}等于当前临时系数{Bn};
所述初始化步骤还包括:设置温度初值A,设置降温速率α,设置温度最小值A0
所述检验步骤中,如果不满足第一条件,进入所述干扰步骤之前,还包括:判断接受当前解的概率是否大于随机数p,如果否,所述当前临时系数{Bn}等于前一临时系数{Bn″},其中E(当前临时系数)-E(前一临时系数)具体为所述第一类方程的空间偏导数在利用当前临时系数{Bn}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果与所述第一类方程的空间偏导数在利用前一临时系数{Bn″}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果之差,所述随机数p具体为0到1之间的随机数;
所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:判断所述A是否大于A0
如果A大于A0,则A=A*α,重新设置优化系数输出条件,重新进入所述干扰步骤;
如果A小于等于A0,进入所述输出步骤。
3.根据权利要求2所述的方法,其特征在于,
所述计算步骤还包括:设置当前离散值为无解状态;
所述获取步骤还包括:设置当前离散值为有解状态,判断是否所述当前离散值<π,如果是,将所述当前离散值增加一个离散间隔作为当前离散值,重新进入所述计算步骤,如果否,进入所述输出步骤;
所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:如果所述A小于等于A0,则判断是否所述当前离散值<π,如果所述当前离散值<π,并且当前离散值为有解状态,则将所述当前离散值增加一个离散间隔作为当前离散值,重新进入所述计算步骤;
如果所述当前离散值≥π或者当前离散值为无解状态,进入所述输出步骤。
4.根据权利要求1所述的方法,其特征在于,还包括:
通过计算所述第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C与所述第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的结果之差,获得第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的误差。
5.根据权利要求4所述的方法,其特征在于,所述将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn},具体为将第一类待选结果中精度覆盖范围最大,且误差和最小的当前临时系数{Bn}作为第一类优化系数{bn};
所述当前临时系数{Bn}的误差和通过计算所述第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的误差之和获得。
6.根据权利要求1所述的方法,其特征在于,还包括:
当所述第一类方程为一阶偏微分方程,所述有限差分格式不是交错网格有限差分时,限定第一类优化系数{bn}满足优化条件,所述优化条件包括:
限定所述当前临时系数{Bn}包括第一类临时系数{B-m}、中间临时系数B0和第二类临时系数{Bm},其中m>0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}以中间临时系数B0为中心奇对称;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
限定所述当前临时系数{Bn}的总和为0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,越邻近中间临时系数B0的系数的绝对值越大;
当所述第一类方程为二阶偏微分方程,所述有限差分格式不是交错网格有限差分时,限定第一类优化系数{bn}满足优化条件,所述优化条件包括:
限定所述当前临时系数{Bn}包括第一类临时系数{B-m}、中间临时系数B0和第二类临时系数{Bm},其中m>0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}以中间临时系数B0为中心偶对称;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
限定所述当前临时系数{Bn}的总和为0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,越邻近中间临时系数B0的系数的绝对值越大;
当所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分时,限定第一类优化系数{bn}满足优化条件,所述优化条件包括:
限定所述当前临时系数{Bn}包括第一类临时系数{B-m+1}、中间临时系数B1和第二类临时系数{Bm},其中m>1;
限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}以中间临时系数B1为中心奇对称;
限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}中,越邻近中间临时系数B1的系数的绝对值越大。
7.根据权利要求6所述的方法,其特征在于,
当所述第一类方程为一阶或二阶偏微分方程,所述有限差分格式不是交错网格有限差分时,所述计算步骤的随机产生至少一组当前临时系数{Bn},具体通过以下步骤产生:
对应每个待求的第二类临时系数Bm各分配一个第一类随机数rm,其中0≤rm≤1;
根据计算出所述第二类临时系数{Bm}的值,其中为Bm预设的浮动上限,为Bm预设的浮动下限;
根据第一类优化系数{bn}的优化条件和所述第二类临时系数{Bm}的值,求出第一类临时系数{B-m}与中间临时系数B0的值;
当所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分时,所述计算步骤的随机产生至少一组当前临时系数{Bn},具体通过以下步骤产生:
对应每个待求的第二类临时系数Bm各分配一个第一类随机数rm,其中0≤rm≤1;
根据计算出所述第二类临时系数{Bm}的值,其中为Bm预设的浮动上限,为Bm预设的浮动下限;
根据第一类优化系数{bn}的优化条件和所述第二类临时系数{Bm}的值,求出第一类临时系数{B-m+1}与中间临时系数B1的值。
8.根据权利要求1到7任意一项所述的方法,其特征在于,所述优化系数输出条件具体为重新进入干扰步骤的次数超过预设的干扰次数阀值。
9.根据权利要求1所述的方法,其特征在于,所述有限差分格式不是交错网格有限差分时,所述预设的误差限T具体为0.0001。
10.根据权利要求1所述的方法,其特征在于,所述有限差分是交错网格有限差分时,所述预设的误差限T具体为0.00005。
11.根据权利要求1所述的方法,其特征在于,
当所述第一类方程为一阶偏微分方程,所述有限差分格式不是交错网格有限差分时,所述判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( K x ( i ) , T ) ≡ m a x 0 ≤ k x ( i ) | - K x ( i ) Δ - Σ n = - N / 2 N / 2 B n s i n ( - K x ( i ) Δ n ) | ≤ T ;
当所述第一类方程为二阶偏微分方程,所述有限差分格式不是交错网格有限差分时,所述判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( K x ( i ) , T ) ≡ m a x 0 ≤ k x ( i ) | - K x ( i ) 2 Δ 2 - Σ n = - N / 2 N / 2 B n c o s ( nK x ( i ) Δ ) | ≤ T ;
当所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分时,所述判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( K x ( i ) , T ) ≡ m a x 0 ≤ k x ( i ) | - K x ( i ) Δ - Σ n = - N / 2 N / 2 b n sin [ ( 0.5 - n ) K x ( i ) Δ ] | ≤ T ;
其中Δ为震源点速度模型的空间网格间距。
12.根据权利要求1所述的方法,其特征在于,
当所述第一类方程具体为一阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中0.0834≤b-2≤0.1985,-0.1985≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中-0.0357≤b-3≤-0.0167,0.1501≤b-2≤0.2912,-0.2912≤b2≤-0.1501,0.0167≤b3≤0.0357;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中0.0036≤b-4≤0.0097,-0.0669≤b-3≤-0.0381,0.2001≤b-2≤0.3698,-0.3698≤b2≤-0.2001,0.0381≤b3≤0.0669,-0.0097≤b4≤-0.0036;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中-0.0078≤b-5≤-0.0008,0.01≤b-4≤0.0299,-0.1337≤b-3≤-0.0596,0.2381≤b-2≤0.3325,-0.3325≤b2≤-0.2381,0.0596≤b3≤0.1337,-0.0299≤b4≤-0.01,0.0008≤b5≤0.0078;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中0.0001≤b-6≤0.0071,-0.0148≤b-5≤-0.0026,0.0179≤b-4≤0.0588,-0.1527≤b-3≤-0.0794,0.2679≤b-2≤0.3766,-0.3766≤b2≤-0.2679,0.0794≤b3≤0.1527,-0.0588≤b4≤-0.0179,0.0026≤b5≤0.0148,-0.0071≤b6≤-0.0001;
当所述第一类方程具体为一阶偏微分方程,所述有限差分是交错网格有限差分时,
用于控制4阶交错网格的有限差分格式的第一类优化系数bn具体为b-1,b1,b2,其中0.04167≤b-1≤0.0913,-0.0913≤b2≤-0.04167;
用于控制6阶交错网格的有限差分格式的第一类优化系数bn具体为b-2,b-1,b1,b2,b3,其中-0.0761≤b-2≤-0.0047,0.0652≤b-1≤0.1820,-0.1820≤b2≤-0.0652,0.0047≤b3≤0.0761;
用于控制8阶交错网格的有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b1,b2,b3,b4,其中0.0007≤b-3≤0.0034,-0.0188≤b-2≤-0.0096,0.0798≤b-1≤0.1465,-0.1465≤b2≤-0.0798,0.0096≤b3≤0.0188,-0.0034≤b4≤-0.0007;
用于控制10阶交错网格的有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,其中-0.0088≤b-4≤-0.0002,0.0018≤b-3≤0.0084,-0.0139≤b-2≤-0.0298,0.0898≤b-1≤0.1969,-0.1969≤b2≤-0.0898,0.0139≤b3≤0.0298,-0.0084≤b4≤-0.0018,0.0002≤b5≤0.0088;
用于控制12阶交错网格的有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,b6,其中0.0002≤b-5≤0.009,-0.0046≤b-4≤-0.0004,0.0030≤b-3≤0.0979,-0.0599≤b-2≤-0.0175,0.0970≤b-1≤0.1953,-0.1953≤b2≤-0.0970,0.0175≤b3≤0.0599,-0.0979≤b4≤-0.0030,0.0004≤b5≤0.0046,-0.009≤b6≤-0.0002;
当所述第一类方程具体为二阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中-0.1648≤b-2≤-0.0834,-0.1648≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中0.0112≤b-3≤0.0373,-0.3018≤b-2≤-0.1510,-0.3018≤b2≤-0.1510,0.0112≤b3≤0.0373;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中-0.0086≤b-4≤-0.0018,0.0254≤b-3≤0.0585,-0.3855≤b-2≤-0.2001,-0.3855≤b2≤-0.2001,0.0254≤b3≤0.0585,-0.0086≤b4≤-0.0018;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中0.0004≤b-5≤0.0038,-0.0188≤b-4≤-0.0050,0.0397≤b-3≤0.0837,-0.4826≤b-2≤-0.2384,-0.4826≤b2≤-0.2384,0.0397≤b3≤0.0837,-0.0188≤b4≤-0.0050,0.0004≤b5≤0.0038;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中-0.0037≤b-6≤-0.0007,0.0011≤b-5≤0.0077,-0.0327≤b-4≤-0.0090,0.0530≤b-3≤0.1128,-0.3927≤b-2≤-0.2679,-0.3927≤b2≤-0.2679,0.0530≤b3≤0.1128,-0.0327≤b4≤-0.0090,0.0011≤b5≤0.0077,-0.0037≤b6≤-0.0007。
13.一种基于优化系数的地震波场模拟装置,其特征在于,包括:
预处理单元:用于获取震源点激发的波动数据,所述震源点激发的波动数据至少包括模型介质的波动速度、震源点空间坐标和震源点时间坐标;获取震源点激发的地震波场模拟涉及的第一类方程;
模拟单元:用于将所述震源点激发的波动数据作为所述第一类方程的输入数据,应用一种优化系数获取方法获取的第一类优化系数{bn}控制有限差分格式对震源点激发的地震波场进行模拟;
所述优化系数获取方法包括初始化步骤、计算步骤、检验步骤、获取步骤、干扰步骤和输出步骤:
所述初始化步骤包括:
设置误差限T的值;
设置当前离散值的初值;
设置优化系数输出条件;
所述计算步骤包括:
随机产生至少一组当前临时系数{Bn},其中 为Bn预设的浮动上限,为Bn预设的浮动下限,其中所述当前临时系数{Bn}中Bn的个数由有限差分格式具体采用的阶数N决定;
所述检验步骤包括:
判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件;
其中,所述第一条件具体为理想值与实际值之间的差值E小于或者等于预设的误差限T,所述理想值具体为第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C,所述实际值具体为所述第一类方程的空间偏导数在利用当前临时系数Bn控制的有限差分格式的傅里叶变换在离散变量Kx(i)取第i个离散值时的结果,所述离散变量Kx(i)的离散值的范围为0≤Kx(i)<π,C为所述第一类方程的空间偏导数的阶数,为虚数单位;
如果满足第一条件,进入所述获取步骤;
如果不满足第一条件,进入所述干扰步骤;
所述获取步骤包括:
将所述当前临时系数{Bn}加入第一类待选结果;
根据判断所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值是否均满足第一条件,获取所述当前临时系数{Bn}的精度覆盖范围,所述精度覆盖范围具体为所述当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)取所述精度覆盖范围内的任意离散值均满足第一条件的最大离散值;
所述干扰步骤包括:
判断优化系数输出条件是否满足;
如果优化系数输出条件未满足,将所述当前临时系数{Bn}在当前基础上进行调整,所述当前临时系数{Bn}调整后的值不超过{Bn}预设的浮动上限和下限,更新所述当前临时系数{Bn}为当前临时系数{Bn}调整后的值,进入所述检验步骤;
如果优化系数输出条件满足,进入所述输出步骤;
所述输出步骤包括:
将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn}。
14.根据权利要求13所述的装置,其特征在于,
在所述计算步骤中随机产生一组当前临时系数{Bn};
在所述计算步骤之后,进入检验步骤之前,还包括:将当前临时系数{Bn}的值在当前基础上进行调整,调整后的值不超过{Bn}预设的浮动上限和下限,获得调整后临时系数{Bn′};
前一临时系数{Bn″}等于当前临时系数{Bn};
所述当前临时系数{Bn}等于所述调整后临时系数{Bn′};
在所述获取步骤中,还包括:所述前一临时系数{Bn″}等于当前临时系数{Bn};
所述初始化步骤还包括:设置温度初值A,设置降温速率α,设置温度最小值A0
所述检验步骤中,如果不满足第一条件,进入所述干扰步骤之前,还包括:判断接受当前解的概率是否大于随机数p,如果否,所述当前临时系数{Bn}等于前一临时系数{Bn″},其中E(当前临时系数)-E(前一临时系数)具体为所述第一类方程的空间偏导数在利用当前临时系数{Bn}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果与所述第一类方程的空间偏导数在利用前一临时系数{Bn″}控制的有限差分格式的傅立叶变换在离散变量取当前离散值时的结果之差,所述随机数p具体为0到1之间的随机数;
所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:判断所述A是否大于A0
如果A大于A0,则A=A*α,重新设置优化系数输出条件,重新进入所述干扰步骤;
如果A小于等于A0,进入所述输出步骤。
15.根据权利要求14所述的装置,其特征在于,
所述计算步骤还包括:设置当前离散值为无解状态;
所述获取步骤还包括:设置当前离散值为有解状态,判断是否所述当前离散值<π,如果是,将所述当前离散值增加一个离散间隔作为当前离散值,重新进入所述计算步骤,如果否,进入所述输出步骤;
所述干扰步骤中,如果优化系数输出条件满足,进入所述输出步骤之前还包括:如果所述A小于等于A0,则判断是否所述当前离散值<π,如果所述当前离散值<π,并且当前离散值为有解状态,则将所述当前离散值增加一个离散间隔作为当前离散值,重新进入所述计算步骤;
如果所述当前离散值≥π或者当前离散值为无解状态,进入所述输出步骤。
16.根据权利要求13所述的装置,其特征在于,还包括:
通过计算所述第一类方程的空间偏导数的傅里叶变换的结果(jKx(i))C与所述第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的结果之差,获得第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的误差。
17.根据权利要求16所述的装置,其特征在于,所述将第一类待选结果中精度覆盖范围最大的当前临时系数{Bn}作为第一类优化系数{bn},具体为将第一类待选结果中精度覆盖范围最大,且误差和最小的当前临时系数{Bn}作为第一类优化系数{bn};
所述当前临时系数{Bn}的误差和通过计算所述第一类待选结果中每个当前临时系数{Bn}控制的有限差分格式的傅里叶变换在离散变量Kx(i)取精度覆盖范围内每个离散值时的误差之和获得。
18.根据权利要求13所述的装置,其特征在于,还包括:
当所述第一类方程为一阶偏微分方程,所述有限差分格式不是交错网格有限差分时,限定第一类优化系数{bn}满足优化条件,所述优化条件包括:
限定所述当前临时系数{Bn}包括第一类临时系数{B-m}、中间临时系数B0和第二类临时系数{Bm},其中m>0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}以中间临时系数B0为中心奇对称;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
限定所述当前临时系数{Bn}的总和为0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,越邻近中间临时系数B0的系数的绝对值越大;
当所述第一类方程为二阶偏微分方程,所述有限差分格式不是交错网格有限差分时,限定第一类优化系数{bn}满足优化条件,所述优化条件包括:
限定所述当前临时系数{Bn}包括第一类临时系数{B-m}、中间临时系数B0和第二类临时系数{Bm},其中m>0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}以中间临时系数B0为中心偶对称;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
限定所述当前临时系数{Bn}的总和为0;
限定所述第一类临时系数{B-m}和第二类临时系数{Bm}中,越邻近中间临时系数B0的系数的绝对值越大;
当所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分时,限定第一类优化系数{bn}满足优化条件,所述优化条件包括:
限定所述当前临时系数{Bn}包括第一类临时系数{B-m+1}、中间临时系数B1和第二类临时系数{Bm},其中m>1;
限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}以中间临时系数B1为中心奇对称;
限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}中,相邻系数相乘结果为负数;
限定所述第一类临时系数{B-m+1}和第二类临时系数{Bm}中,越邻近中间临时系数B1的系数的绝对值越大。
19.根据权利要求18所述的装置,其特征在于,
当所述第一类方程为一阶或二阶偏微分方程,所述有限差分格式不是交错网格有限差分时,所述计算步骤的随机产生至少一组当前临时系数{Bn},具体通过以下步骤产生:
对应每个待求的第二类临时系数Bm各分配一个第一类随机数rm,其中0≤rm≤1;
根据计算出所述第二类临时系数{Bm}的值,其中为Bm预设的浮动上限,为Bm预设的浮动下限;
根据第一类优化系数{bn}的优化条件和所述第二类临时系数{Bm}的值,求出第一类临时系数{B-m}与中间临时系数B0的值;
当所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分时,所述计算步骤的随机产生至少一组当前临时系数{Bn},具体通过以下步骤产生:
对应每个待求的第二类临时系数Bm各分配一个第一类随机数rm,其中0≤rm≤1;
根据计算出所述第二类临时系数{Bm}的值,其中为Bm预设的浮动上限,为Bm预设的浮动下限;
根据第一类优化系数{bn}的优化条件和所述第二类临时系数{Bm}的值,求出第一类临时系数{B-m+1}与中间临时系数B1的值。
20.根据权利要求13到19任意一项所述的装置,其特征在于,所述优化系数输出条件具体为重新进入干扰步骤的次数超过预设的干扰次数阀值。
21.根据权利要求13所述的装置,其特征在于,所述有限差分格式不是交错网格有限差分时,所述预设的误差限T具体为0.0001。
22.根据权利要求13所述的装置,其特征在于,所述有限差分是交错网格有限差分时,所述预设的误差限T具体为0.00005。
23.根据权利要求13所述的装置,其特征在于,
当所述第一类方程为一阶偏微分方程,所述有限差分格式不是交错网格有限差分时,所述判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( K x ( i ) , T ) ≡ m a x 0 ≤ k x ( i ) | - K x ( i ) Δ - Σ n = - N / 2 N / 2 B n s i n ( - K x ( i ) Δ n ) | ≤ T ;
当所述第一类方程为二阶偏微分方程,所述有限差分格式不是交错网格有限差分时,所述判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( K x ( i ) , T ) ≡ m a x 0 ≤ k x ( i ) | - K x ( i ) 2 Δ 2 - Σ n = - N / 2 N / 2 B n c o s ( nK x ( i ) Δ ) | ≤ T ;
当所述第一类方程为一阶偏微分方程,所述有限差分格式是交错网格有限差分时,所述判断当前临时系数{Bn}控制的有限差分格式的离散变量Kx(i)从0到当前离散值,是否均满足第一条件,具体利用以下目标函数进行判断:
E ( K x ( i ) , T ) ≡ m a x 0 ≤ k x ( i ) | - K x ( i ) Δ - Σ n = - N / 2 N / 2 b n sin [ ( 0.5 - n ) K x ( i ) Δ ] | ≤ T ;
其中Δ为震源点速度模型的空间网格间距。
24.根据权利要求13所述的装置,其特征在于,
当所述第一类方程具体为一阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中0.0834≤b-2≤0.1985,-0.1985≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中-0.0357≤b-3≤-0.0167,0.1501≤b-2≤0.2912,-0.2912≤b2≤-0.1501,0.0167≤b3≤0.0357;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中0.0036≤b-4≤0.0097,-0.0669≤b-3≤-0.0381,0.2001≤b-2≤0.3698,-0.3698≤b2≤-0.2001,0.0381≤b3≤0.0669,-0.0097≤b4≤-0.0036;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中-0.0078≤b-5≤-0.0008,0.01≤b-4≤0.0299,-0.1337≤b-3≤-0.0596,0.2381≤b-2≤0.3325,-0.3325≤b2≤-0.2381,0.0596≤b3≤0.1337,-0.0299≤b4≤-0.01,0.0008≤b5≤0.0078;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中0.0001≤b-6≤0.0071,-0.0148≤b-5≤-0.0026,0.0179≤b-4≤0.0588,-0.1527≤b-3≤-0.0794,0.2679≤b-2≤0.3766,-0.3766≤b2≤-0.2679,0.0794≤b3≤0.1527,-0.0588≤b4≤-0.0179,0.0026≤b5≤0.0148,-0.0071≤b6≤-0.0001;
当所述第一类方程具体为一阶偏微分方程,所述有限差分是交错网格有限差分时,
用于控制4阶交错网格的有限差分格式的第一类优化系数bn具体为b-1,b1,b2,其中0.04167≤b-1≤0.0913,-0.0913≤b2≤-0.04167;
用于控制6阶交错网格的有限差分格式的第一类优化系数bn具体为b-2,b-1,b1,b2,b3,其中-0.0761≤b-2≤-0.0047,0.0652≤b-1≤0.1820,-0.1820≤b2≤-0.0652,0.0047≤b3≤0.0761;
用于控制8阶交错网格的有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b1,b2,b3,b4,其中0.0007≤b-3≤0.0034,-0.0188≤b-2≤-0.0096,0.0798≤b-1≤0.1465,-0.1465≤b2≤-0.0798,0.0096≤b3≤0.0188,-0.0034≤b4≤-0.0007;
用于控制10阶交错网格的有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,其中-0.0088≤b-4≤-0.0002,0.0018≤b-3≤0.0084,-0.0139≤b-2≤-0.0298,0.0898≤b-1≤0.1969,-0.1969≤b2≤-0.0898,0.0139≤b3≤0.0298,-0.0084≤b4≤-0.0018,0.0002≤b5≤0.0088;
用于控制12阶交错网格的有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b1,b2,b3,b4,b5,b6,其中0.0002≤b-5≤0.009,-0.0046≤b-4≤-0.0004,0.0030≤b-3≤0.0979,-0.0599≤b-2≤-0.0175,0.0970≤b-1≤0.1953,-0.1953≤b2≤-0.0970,0.0175≤b3≤0.0599,-0.0979≤b4≤-0.0030,0.0004≤b5≤0.0046,-0.009≤b6≤-0.0002;
当所述第一类方程具体为二阶偏微分方程,所述有限差分不是交错网格有限差分时,
用于控制4阶有限差分格式的第一类优化系数bn具体为b-2,b-1,b0,b1,b2,其中-0.1648≤b-2≤-0.0834,-0.1648≤b2≤-0.0834;
用于控制6阶有限差分格式的第一类优化系数bn具体为b-3,b-2,b-1,b0,b1,b2,b3,其中0.0112≤b-3≤0.0373,-0.3018≤b-2≤-0.1510,-0.3018≤b2≤-0.1510,0.0112≤b3≤0.0373;
用于控制8阶有限差分格式的第一类优化系数bn具体为b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,其中-0.0086≤b-4≤-0.0018,0.0254≤b-3≤0.0585,-0.3855≤b-2≤-0.2001,-0.3855≤b2≤-0.2001,0.0254≤b3≤0.0585,-0.0086≤b4≤-0.0018;
用于控制10阶有限差分格式的第一类优化系数bn具体为b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,其中0.0004≤b-5≤0.0038,-0.0188≤b-4≤-0.0050,0.0397≤b-3≤0.0837,-0.4826≤b-2≤-0.2384,-0.4826≤b2≤-0.2384,0.0397≤b3≤0.0837,-0.0188≤b4≤-0.0050,0.0004≤b5≤0.0038;
用于控制12阶有限差分格式的第一类优化系数bn具体为b-6,b-5,b-4,b-3,b-2,b-1,b0,b1,b2,b3,b4,b5,b6,其中-0.0037≤b-6≤-0.0007,0.0011≤b-5≤0.0077,-0.0327≤b-4≤-0.0090,0.0530≤b-3≤0.1128,-0.3927≤b-2≤-0.2679,-0.3927≤b2≤-0.2679,0.0530≤b3≤0.1128,-0.0327≤b4≤-0.0090,0.0011≤b5≤0.0077,-0.0037≤b6≤-0.0007。
CN201210343161.1A 2012-09-14 2012-09-14 一种基于优化系数的地震波场模拟方法及装置 Expired - Fee Related CN103675905B (zh)

Priority Applications (3)

Application Number Priority Date Filing Date Title
CN201210343161.1A CN103675905B (zh) 2012-09-14 2012-09-14 一种基于优化系数的地震波场模拟方法及装置
US14/117,307 US20150134308A1 (en) 2012-09-14 2012-11-05 Method and device for acquiring optimization coefficient, and related method and device for simulating wave field
PCT/CN2012/084083 WO2014040338A1 (zh) 2012-09-14 2012-11-05 一种优化系数获取方法、装置及相关波场模拟方法、装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201210343161.1A CN103675905B (zh) 2012-09-14 2012-09-14 一种基于优化系数的地震波场模拟方法及装置

Publications (2)

Publication Number Publication Date
CN103675905A CN103675905A (zh) 2014-03-26
CN103675905B true CN103675905B (zh) 2016-10-05

Family

ID=50277532

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201210343161.1A Expired - Fee Related CN103675905B (zh) 2012-09-14 2012-09-14 一种基于优化系数的地震波场模拟方法及装置

Country Status (3)

Country Link
US (1) US20150134308A1 (zh)
CN (1) CN103675905B (zh)
WO (1) WO2014040338A1 (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105844000A (zh) * 2016-03-18 2016-08-10 江苏铨铨信息科技有限公司 一种mcc 表面海流反演方法
CN107766085A (zh) * 2017-09-30 2018-03-06 昂纳信息技术(深圳)有限公司 一种提高系统控制稳定性的方法、装置及存储装置
CN110873895A (zh) * 2018-08-31 2020-03-10 中国石油化工股份有限公司 一种变网格微地震逆时干涉定位方法
CN109270575B (zh) * 2018-11-02 2019-11-26 河南理工大学 一种基于建筑物地震响应等效的爆破地震波模型构造方法
CN112907885B (zh) * 2021-01-12 2022-08-16 中国计量大学 基于scnn的分散集中式家居图像火灾报警系统及方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101013161A (zh) * 2007-01-15 2007-08-08 中国石油大港油田勘探开发研究院 基于叠前波场模拟的地震勘探层位标定方法
CN101576621A (zh) * 2008-05-07 2009-11-11 王振华 海底电缆双检地震勘探的数据处理方法及数据处理装置

Family Cites Families (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB2329043B (en) * 1997-09-05 2000-04-26 Geco As Method of determining the response caused by model alterations in seismic simulations
GB2381314B (en) * 2001-10-26 2005-05-04 Westerngeco Ltd A method of and an apparatus for processing seismic data
FR2843202B1 (fr) * 2002-08-05 2004-09-10 Inst Francais Du Petrole Methode pour former un modele representatif de la distribution d'une grandeur physique dans une zone souterraine, affranchi de l'effet de bruits correles entachant des donnees d'exploration
US8457899B2 (en) * 2007-12-14 2013-06-04 Shell Oil Company Method of processing data obtained from seismic prospecting
JP5279016B2 (ja) * 2008-11-21 2013-09-04 国立大学法人 東京大学 地球内部を伝播する地震波についての数値解析に用いる演算子生成方法、及び演算子生成装置、並びに地球内部を伝播する地震波についての数値解析を行うシミュレーション装置
WO2015042386A1 (en) * 2013-09-20 2015-03-26 Westerngeco Llc Eikonal solver for quasi p-waves in anisotropic media
US10185046B2 (en) * 2014-06-09 2019-01-22 Exxonmobil Upstream Research Company Method for temporal dispersion correction for seismic simulation, RTM and FWI
US20160320509A1 (en) * 2015-04-30 2016-11-03 Saudi Arabian Oil Company Suppressing near-surface scattered surface waves

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101013161A (zh) * 2007-01-15 2007-08-08 中国石油大港油田勘探开发研究院 基于叠前波场模拟的地震勘探层位标定方法
CN101576621A (zh) * 2008-05-07 2009-11-11 王振华 海底电缆双检地震勘探的数据处理方法及数据处理装置

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
从瞬变电磁场到波场的优化算法;李貅;《地球物理学报》;20050930;第48卷(第5期);1185-1190 *
基于多参量的模拟退火全局优化傅里叶有限差分算子;朱遂伟;《地球物理学报》;20081130;第51卷(第6期);1844-1850 *
高阶优化傅里叶有限差分算子偏移;朱遂伟;《石油地球物理勘探》;20091231;第44卷(第6期);680-684 *

Also Published As

Publication number Publication date
WO2014040338A1 (zh) 2014-03-20
CN103675905A (zh) 2014-03-26
US20150134308A1 (en) 2015-05-14

Similar Documents

Publication Publication Date Title
Sjögreen et al. A fourth order accurate finite difference scheme for the elastic wave equation in second order formulation
Nilsson et al. Stable difference approximations for the elastic wave equation in second order formulation
CN106842306B (zh) 一种全局优化的交错网格有限差分正演模拟方法和装置
Liu et al. Prediction of broadband ground-motion time histories: Hybrid low/high-frequency method with correlated random source parameters
CN103675905B (zh) 一种基于优化系数的地震波场模拟方法及装置
CN102057368B (zh) 利用最大连续场在三维体积模型中分布性质
CN107066772B (zh) 非平稳地震作用下桥梁系统碰撞间隙宽度的概率评价方法
Minisini et al. Local time stepping with the discontinuous Galerkin method for wave propagation in 3D heterogeneous media
Chung et al. A staggered discontinuous Galerkin method for the simulation of seismic waves with surface topography
CN108984939A (zh) 基于3D Gauss-FFT的三维重力场正演方法
CN110598367A (zh) 一种足迹引导的高效航空电磁法数值模拟方法
Song Developing a generalized pseudo-dynamic source model of M w 6.5–7.0 to simulate strong ground motions
Alkhalifah et al. An eikonal-based formulation for traveltime perturbation with respect to the source location
Gadylshin et al. Optimization of the training dataset for numerical dispersion mitigation neural network
Ba et al. Near-fault broadband seismograms synthesis in a stratified transversely isotropic half-space using a semi-analytical frequency-wavenumber method
Nakai et al. Observation site selection for physical model parameter estimation towards process-driven seismic wavefield reconstruction
CN105676280A (zh) 基于旋转交错网格的双相介质地质数据获取方法和装置
CN109239776A (zh) 一种地震波传播正演模拟方法和装置
Morikawa et al. Shaking maps for scenario earthquakes by applying the upgraded version of the strong ground motion prediction method “recipe”
Restrepo et al. Virtual topography: A fictitious domain approach for analyzing free‐surface irregularities in large‐scale earthquake ground motion simulation
Magrin Multi-scale seismic hazard scenarios
Jing et al. Optimization algorithm for rapid 3D gravity inversion
Gallovic et al. High-frequency directivity in strong ground motion modeling methods
MX2011003852A (es) Atributos de procesamiento de imagen con inversion de tiempo.
Ulrich et al. Rapidness and robustness of finite-source inversion of the 2011 M w 9.0 Tohoku earthquake by an elliptical-patches method using continuous GPS and acceleration data

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20161005

Termination date: 20190914