一种基于随机梯度采样的弹性波全波形反演方法技术

本发明专利技术涉及一种基于随机梯度采样的弹性波全波形反演方法，包括以下步骤：步骤1)对原始地震数据进行反演前预处理；步骤2)对预处理后的地震数据进行傅里叶变换，得到频率域真实地震记录；步骤3)建立弹性参数模型并初始化，设定弹性波全波形的反演参数并定义观测系统；步骤4)采用随机梯度采样弹性波全波形反演方法，基于频率域真实地震记录进行弹性参数模型更新，得到最终的弹性参数模型，完成反演。与现有技术相比，本发明专利技术能够实现在较差的初始弹性参数模型或缺乏低频地震数据下有效进行弹性波全波形反演，一定程度上减弱周期跳跃现象对反演结果带来的不利影响，具有实现算法简单、不额外增加内存消耗与计算量等优点。不额外增加内存消耗与计算量等优点。不额外增加内存消耗与计算量等优点。

【技术实现步骤摘要】
一种基于随机梯度采样的弹性波全波形反演方法


[0001]本专利技术涉及勘探地震学中的速度建模领域，尤其是涉及一种基于随机梯度采样的弹性波全波形反演方法。

技术介绍

[0002]弹性波全波形反演可以利用多分量弹性波地震数据来获取地下弹性参数。随着近年来计算机水平与地震采集技术的高速发展，一定程度上解决了弹性波全波形反演大幅的计算量与内存消耗问题，成功应用于大量实际数据。然而，弹性波全波形反演在获得了丰富的地下物性参数信息的同时，也需要迎接严峻的挑战。由于全波形反演是一个高度非线性且非凸优化问题，使得反演容易陷入局部极值。较低的横波速度会加剧弹性波全波形反演问题的非线性程度，相比纵波速度，反演横波速度需要构建更加准确的初始模型来避免陷入局部极值。同时，低频地震数据的缺失对重建横波速度模型更加不利，“周期跳跃”现象会更加凸显。另一方面，由于弹性波在地下介质传播机制复杂，会产生绕射、转换波等现象，参数耦合效应也会给反演结果带来十分不利的影响。
[0003]针对弹性波全波形反演的参数耦合效应，大量学者提出了一系列的解决方案。诸如多步反演策略来逐步恢复纵波速度(阻抗)、横波速度(阻抗)、密度等，参数化方式的合理选择，P/S波场分离，可以有效地消除纵波速度扰动在反演结果上的混叠干扰，实现参数解耦。Xu和McMechan提出了一种多步长的反演方法，来约束强弱不同参数的更新，也能部分地压制参数耦合效应。
[0004]然而，对于“周期跳跃”问题的研究还较少。多尺度反演策略(从低频数据到高频数据)虽然能优先获得背景信息，避免高波数信息的过早恢复，但是原则上没有解决低频数据缺失问题。其他如利用走时层析提供较好的背景初始模型，或使用地震波形包络信息、相位信息反演等，利用地下散射角滤波的方式获取参数背景、最佳输运距离弹性波反演等方法也能一定程度上处理“周期跳跃”问题，改善反演结果并获得较高精度与分辨率的参数模型。
[0005]但是，上述这些方法均改变了原始弹性波全波形反演的目标泛函，增加了算法复杂度，且有些算法在实施过程中涉及大量额外的计算和内存消耗，阻碍了其实际应用。

技术实现思路

[0006]本专利技术的目的就是为了提供一种基于随机梯度采样的弹性波全波形反演方法，实现在较差的初始弹性参数模型或缺乏低频地震数据下有效进行弹性波全波形反演，减弱周期跳跃现象对反演结果带来的不利影响，并简化算法复杂度，提高计算速度。
[0007]本专利技术的目的可以通过以下技术方案来实现：
[0008]一种基于随机梯度采样的弹性波全波形反演方法，包括以下步骤：
[0009]步骤1)对原始地震数据进行反演前预处理；
[0010]步骤2)对预处理后的地震数据进行傅里叶变换，得到频率域真实地震记录；
[0011]步骤3)建立弹性参数模型并初始化，设定弹性波全波形的反演参数并定义观测系统；
[0012]步骤4)采用随机梯度采样弹性波全波形反演方法，基于频率域真实地震记录进行弹性参数模型更新，得到最终的弹性参数模型，完成反演。
[0013]所述弹性参数包括纵波速度、横波速度与密度。
[0014]所述反演参数包括反演频率。
[0015]所述步骤4)包括以下步骤：
[0016]步骤41)执行频率组循环；
[0017]步骤42)在每一个频率组循环内部，执行迭代循环，在上一次迭代更新的弹性参数模型下对每一炮执行炮点端弹性波正传波场模拟，对模拟得到的地表地震记录与波场分别进行傅里叶变换，得到频率域模拟地震记录与模拟波场{u(x,z,ω)＝(u
x
,u
z
)
T
}
n
,(n＝1,2,3
…
Ns)，其中，g表示检波器的空间位置，ω为角频率，Ns为总炮数；
[0018]步骤43)由频率域真实地震记录d
obs
(g,ω)与模拟地震记录d
cal
(g,ω)计算残差地震记录res(g,ω)＝(res
x
,res
z
)
T
，并计算L2范数下真实地震记录与模拟地震记录的残差作为目标函数，其中，为步骤2)获得的频率域真实地震记录；
[0019]步骤44)基于目标函数判断反演是否满足收敛或停止迭代条件，若满足，则跳出当前频率组循环，返回步骤41)，进入下一个频率组循环，直到获得最终的反演结果；若不满足，则执行步骤45)；
[0020]步骤45)对每一炮，将残差记录作为弹性波波场模拟的伴随震源，执行检波端弹性波反传波场模拟，对模拟得到的地表地震记录进行傅里叶变换，得到频率域模拟反传地震记录
[0021]步骤46)采用基于随机梯度采样的计算方法计算目标函数对各弹性参数的梯度；
[0022]步骤47)利用共轭梯度法，进行弹性参数模型更新；
[0023]步骤48)返回步骤42)，进入下一次迭代。
[0024]所述步骤42)具体包括以下步骤：
[0025]步骤421)利用空间四阶、时间二阶有限差分方法，在边界处实施CPML吸收边界条件，对每一炮执行炮点端弹性波正传波场模拟，满足各向同性介质下弹性波位移应力方程：
[0026][0027]其中，u
x
与u
z
为质点位移的X分量与Z分量；σ
xx
，σ
xz
，σ
zx
，σ
zz
为应力分量；λ,μ为拉梅系数，λ＝(α2‑
2β2)ρ,μ＝β2ρ，α,β,ρ分别是纵波速度，横波速度与密度；f
x
与f
z
为体波震源的X分量与Z分量；
[0028]步骤422)对模拟得到的地表地震记录与波场分别进行傅里叶变换，得到频率域模拟地震记录{d
cal
(g,ω)}
n
与模拟波场{u(x,z,ω)}
n
,(n＝1,2,3
…
Ns)。
[0029]所述步骤43)具体包括以下步骤：
[0030]步骤431)计算残差地震记录：
[0031]res(g,ω)＝d
cal
(g,ω)
‑
d
obs
(g,ω)
[0032]其中，d
cal
(g,ω),d
obs
(g,ω)分别为第n炮在空间位置g接收的频率域模拟地震记录与真实地震记录；
[0033]步骤432)计算L2范数下真实地震记录与模拟地震记录的残差作为目标函数：
[0034][0035]其中，res
x,n
(g,ω),res
z,n
(g,ω)分别为第n炮在空间位置g的残差地震记录X分量与Z分量；Nr为第n炮的总检波点数；为共轭。
[0036]所述步骤45)具体包括以下步骤：
[0037]步骤451)利用空本文档来自技高网...

【技术保护点】

【技术特征摘要】
1.一种基于随机梯度采样的弹性波全波形反演方法，其特征在于，包括以下步骤：步骤1)对原始地震数据进行反演前预处理；步骤2)对预处理后的地震数据进行傅里叶变换，得到频率域真实地震记录；步骤3)建立弹性参数模型并初始化，设定弹性波全波形的反演参数并定义观测系统；步骤4)采用随机梯度采样弹性波全波形反演方法，基于频率域真实地震记录进行弹性参数模型更新，得到最终的弹性参数模型，完成反演。2.根据权利要求1所述的一种基于随机梯度采样的弹性波全波形反演方法，其特征在于，所述弹性参数包括纵波速度、横波速度与密度。3.根据权利要求1所述的一种基于随机梯度采样的弹性波全波形反演方法，其特征在于，所述反演参数包括反演频率。4.根据权利要求1所述的一种基于随机梯度采样的弹性波全波形反演方法，其特征在于，所述步骤4)包括以下步骤：步骤41)执行频率组循环；步骤42)在每一个频率组循环内部，执行迭代循环，在上一次迭代更新的弹性参数模型下对每一炮执行炮点端弹性波正传波场模拟，对模拟得到的地表地震记录与波场分别进行傅里叶变换，得到频率域模拟地震记录与模拟波场{u(x,z,ω)＝(u
x
,u
z
)
T
}
n
,(n＝1,2,3
…
Ns)，其中，g表示检波器的空间位置，ω为角频率，Ns为总炮数；步骤43)由频率域真实地震记录d
obs
(g,ω)与模拟地震记录d
cal
(g,ω)计算残差地震记录res(g,ω)＝(res
x
,res
z
)
T
，并计算L2范数下真实地震记录与模拟地震记录的残差作为目标函数，其中，为步骤2)获得的频率域真实地震记录；步骤44)基于目标函数判断反演是否满足收敛或停止迭代条件，若满足，则跳出当前频率组循环，返回步骤41)，进入下一个频率组循环，直到获得最终的反演结果；若不满足，则执行步骤45)；步骤45)对每一炮，将残差记录作为弹性波波场模拟的伴随震源，执行检波端弹性波反传波场模拟，对模拟得到的地表地震记录进行傅里叶变换，得到频率域模拟反传地震记录步骤46)采用基于随机梯度采样的计算方法计算目标函数对各弹性参数的梯度；步骤47)利用共轭梯度法，进行弹性参数模型更新；步骤48)返回步骤42)，进入下一次迭代。5.根据权利要求4所述的一种基于随机梯度采样的弹性波全波形反演方法，其特征在于，所述步骤42)具体包括以下步骤：步骤421)利用空间四阶、时间二阶有限差分方法，在边界处实施CPML吸收边界条件，对每一炮执行炮点端弹性波正传波场模拟，满足各向同性介质下弹性波位移应力方程：
其中，u
x
与u
z
为质点位移的X分量与Z分量；σ
xx
，σ
xz
，σ
zx
，σ
zz
为应力分量；λ,μ为拉梅系数，λ＝(α2‑
2β2)ρ,μ＝β2ρ，α,β,ρ分别是纵波速度，横波速度与密度；f
x
与f
z
为体波震源的X分量与Z分量；步骤422)对模拟得到的地表地震记录与波场分别进行傅里叶变换，得到频率域模拟地震记录{d
cal
(g,ω)}
n
与模拟波场{u(x,z,ω)}
n
,(n＝1,2,3
…
Ns)。6.根据权利要求4所述的一种基于随机梯度采样的弹性波全波形反演方法，其特征在于，所述步骤43)具体包括以下步骤：步骤431)计算残差地震记录：res(g，ω)＝d
cal
(g，ω)
‑
d
obs
(g，ω)其中，d
cal
(g,ω),d
obs
(g,ω)分别为第n炮在空间位置g接收的频率域模拟地震记录与真实地震记录；步骤432)计算L2范数下真实地震记录与模拟地震记录的残差作为目标函数：其中，res
x,n
(g,ω),res
...

【专利技术属性】
技术研发人员：杨积忠，姚宇晖，刘玉柱，
申请(专利权)人：同济大学，
类型：发明
国别省市：