基于最小二乘梯度更新速度模型的地震波全波形反演方法技术

本发明专利技术涉及一种基于最小二乘梯度更新速度模型的地震波全波形反演方法，包括以下步骤:1)获得时间域内地震波场的声波方程；2)构建初始速度模型，设置速度模型更新迭代次数N与允许最小误差值ε；3)构造观测波场数据和计算波场数据的波场误差矢量；4)构造目标函数；5)对目标函数计算获得地震波全波形反演的标准方程；6)引入速度模型更新梯度方向gk和更新量步长α；7)利用最小二乘方法求解速度模型更新梯度方向gk；8)对更新量步长α进行插值计算；9)对速度模型进行更新得：mk＝mk-1+αgk；当|αgk|＜ε时或速度模型更新次数达到速度模型更新迭代次数N时，速度模型更新结束；否则，进入步骤3)。本发明专利技术可快速完成速度模型更新，广泛应用于地震波全波形反演方法中。

【技术实现步骤摘要】

本专利技术涉及一种地震波全波形反演方法，特别是关于一种基于最小二乘梯度更新速度模型的地震波全波形反演方法。
技术介绍
地震波场反演方法很多，如相位反演成像、振幅反演成像以及全波形反演成像。相位反演和振幅反演利用的是波场运动学信息，反演出的速度分布是一种平滑估计值，它是真实速度模型的低频分量，反映的是大尺度的速度结构。基于波动方程理论的全波形反演综合利用了地震记录中振幅、走时和相位等完整的波场信息，通过拟合，全波形反演通常利用给定时窗内的理论波场记录与实测波场记录差值的二次泛函以及额外的地质约束条件建立优化目标函数，利用迭代线性化反演方法或完全非线性反演法求解符合条件的模型，定量提取地下介质的弹性参数，进而为深部大尺度构造演化分析，为勘探地震成像及速度建模等方面提供可靠依据。目前，全波形反演的方式大致可以分为时间域全波形反演与频率域全波形反演。频率域全波形反演理论与时间域全波形反演理论有着良好的对应关系。两种反演方法的观测数据表现形式不同造成正演波场计算方式的差别，而反演思想完全一致，与时间域相比，频率域波场是相互解耦的。不同频率数据对异常体反映能力不一，根据需要既可以利用部分频段数据反演，也可以使用全频段数据同时进行反演。全波形反演可利用逆时偏移产生的成像结果近似来更新速度模型，基于逆时偏移的全波形反演获得低频成像不再需要炮集的低频信息且具有很高的保真度。但是，反演使用大量低频的信息来更新背景速度，而偏移使用高频信息来勾画地下构造的比较精确的边界，由于反演和偏移的目标在不同的频带，因此，在实际处理中，由于噪声和模型假象，高保真度很难达到，且常规全波形反演流程中，每次的全波形反演的非线性迭代过程中，只利用一次的逆时偏移的成像结果作为速度模型更新梯度的方向，这样容易导致更新梯度计算的更新量不足，且收敛速度较慢。
技术实现思路
针对上述问题，本专利技术的目的是提供一种可有效提高速度模型更新梯度的准确性且可快速完成速度模型更新的基于最小二乘梯度更新速度模型的地震波全波形反演方法。为实现上述目的，本专利技术采取以下技术方案：一种基于最小二乘梯度更新速度模型的地震波全波形反演方法，包括以下步骤：1)在时间域内，地震波场的声波方程为：M(x)d2u(x,t)dt2=A(x)u(x,t)+s(x,t);]]>其中，u是波场向量，s是震源向量，x和t分别是空间和时间坐标参数，M为质量矩阵，A为刚度矩阵；2)构建初始速度模型，设置速度模型更新迭代次数N与速度模型允许的最小误差值ε；3)构造观测波场数据和计算波场数据的波场误差矢量Δd；4)利用最小二乘流程和波场误差矢量Δd，构造目标函数C(mk)，目标函数C(mk)的表达式为：C(mk)＝1/2ΔdTΔd；式中，T表示转置；mk为第k次更新的速度模型，更新公式为：mk＝mk-1+Δm；式中，mk-1为第k-1次更新的速度模型；Δm为扰动模型；5)对目标函数进行计算获得地震波全波形反演的标准方程为：Δm=-[∂2C(mk-1)∂mn2]-1∂C(mk-1)∂mn;]]>式中，mn为模型参数；6)引入速度模型更新梯度方向和更新量步长，扰动模型Δm表示为：Δm＝αgk；式中，gk为速度模型更新梯度方向，α为更新量步长；7)利用最小二乘的方法求解地震波全波形反演的速度模型更新梯度方向gk，其公式为：[∂2C(mk-1)∂mn2]gk=-∂C(mk-1)∂mn;]]>8)在速度模型更新梯度方向gk上对更新量步长α进行插值计算求取更新量步长α；9)对速度模型进行更新得：mk＝mk-1+αgk；当|αgk|＜ε时或速度模型更新次数达到速度模型更新迭代次数N时，速度模型更新结束；否则，进入步骤3)。所述步骤3)中，波场误差矢量Δd的计算公式如下：Δd＝dobs-dcal；式中，dobs为地震观测波场数据，dcal为在迭代过程中的计算波场数据；计算波场数据dcal的计算公式如下：式中，为计算波场数据dcal和波场向量u之间的关系函数。所述步骤5)中，地震波全波形反演的标准方程具体计算过程包括以下步骤：(1)对目标函数C(mk)进行二阶泰勒-拉格朗日展开，得到以下形式：C(mk-1+Δm)=C(mk-1)+Σj=1M∂C(mk-1)∂mjΔmj+Σj=1MΣl=1M∂2C(ml-1)∂mj∂mlΔmjΔml+o2(mk-1);]]>式中，j、l分别为二维数据的位置参数；M为正整数；mj、ml分别为更新模型变量mk-1在二维参数模型中某一个方向的分量；(2)对步骤(1)中目标函数C(mk-1+Δm)的展开式以模型参数mn为变量进行求导：∂C(mk)∂mn=∂C(mk-1)∂mn+Σj=1M∂2C(mk-1)∂mj∂mnΔmj;]]>当目标函数C(mk)的导数为零时，目标函数C(mk)得到极值，此时，扰动模型Δm的表达式为：Δm=-[∂2C(mk-1)∂mn2]-1∂C(mk-1)∂mn;]]>扰动模型Δm的表达式即为地震波全波形反演的标准方程。所述步骤8)中，更新量步长α满足伍尔夫条件：C(mk+αgk)<C(mk-1)+αc1gkT▿C(mk);]]>|dT▿C(mk+αgk)/gkT▿C(mk-1)|<c2;]]>式中，c1和c2均为固定参数，其中c1∈(0，1)，c2∈(c1，1)，▽C为目标函数C的梯度。本专利技术由于采取以上技术方案，其具有以下优点：1、本专利技术由于采用迭代非线性反演方法求解符合条件的模型，为深部大尺度构造演化分析、勘探地震成像及速度建模等方面提供可靠依据。2、本专利技术由于采用在传统全波形反演的计算流程基础上，在每次的全波形反演的非线性迭代过程中，引入最小二乘逆时偏移的流程，通过最小二乘的约束和多次本文档来自技高网...

【技术保护点】
一种基于最小二乘梯度更新速度模型的地震波全波形反演方法，包括以下步骤：1)在时间域内，地震波场的声波方程为：M(x)d2u(x,t)dt2=A(x)u(x,t)+s(x,t);]]>其中，u是波场向量，s是震源向量，x和t分别是空间和时间坐标参数，M为质量矩阵，A为刚度矩阵；2)构建初始速度模型，设置速度模型更新迭代次数N与速度模型允许的最小误差值ε；3)构造观测波场数据和计算波场数据的波场误差矢量Δd；4)利用最小二乘流程和波场误差矢量Δd，构造目标函数C(mk)，目标函数C(mk)的表达式为：C(mk)＝1/2ΔdTΔd；式中，T表示转置；mk为第k次更新的速度模型，更新公式为：mk＝mk‑1+Δm；式中，mk‑1为第k‑1次更新的速度模型；Δm为扰动模型；5)对目标函数进行计算获得地震波全波形反演的标准方程为：Δm=-[∂2C(mk-1)∂mn2]-1∂C(mk-1)∂mn;]]>式中，mn为模型参数；6)引入速度模型更新梯度方向和更新量步长，扰动模型Δm表示为：Δm＝αgk；式中，gk为速度模型更新梯度方向，α为更新量步长；7)利用最小二乘的方法求解地震波全波形反演的速度模型更新梯度方向gk，其公式为：[∂2C(mk-1)∂mn2]gk=-∂C(mk-1)∂mn;]]>8)在速度模型更新梯度方向gk上对更新量步长α进行插值计算求取更新量步长α；9)对速度模型进行更新得：mk＝mk‑1+αgk；当|αgk|＜ε时或速度模型更新次数达到速度模型更新迭代次数N时，速度模型更新结束；否则，进入步骤3)。...

【技术特征摘要】
1.一种基于最小二乘梯度更新速度模型的地震波全波形反演方法，包括以下
步骤：
1)在时间域内，地震波场的声波方程为：
M(x)d2u(x,t)dt2=A(x)u(x,t)+s(x,t);]]>其中，u是波场向量，s是震源向量，x和t分别是空间和时间坐标参数，M为
质量矩阵，A为刚度矩阵；
2)构建初始速度模型，设置速度模型更新迭代次数N与速度模型允许的最小
误差值ε；
3)构造观测波场数据和计算波场数据的波场误差矢量Δd；
4)利用最小二乘流程和波场误差矢量Δd，构造目标函数C(mk)，目标函数
C(mk)的表达式为：
C(mk)＝1/2ΔdTΔd；
式中，T表示转置；
mk为第k次更新的速度模型，更新公式为：
mk＝mk-1+Δm；
式中，mk-1为第k-1次更新的速度模型；Δm为扰动模型；
5)对目标函数进行计算获得地震波全波形反演的标准方程为：
Δm=-[∂2C(mk-1)∂mn2]-1∂C(mk-1)∂mn;]]>式中，mn为模型参数；
6)引入速度模型更新梯度方向和更新量步长，扰动模型Δm表示为：
Δm＝αgk；
式中，gk为速度模型更新梯度方向，α为更新量步长；
7)利用最小二乘的方法求解地震波全波形反演的速度模型更新梯度方向gk，
其公式为：
[∂2C(mk-1)∂mn2]gk=-∂C(mk-1)∂mn;]]>8)在速度模型更新梯度方向gk上对更新量步长α进行插值计算求取更新量步
长α；
9)对速度模型进行更新得：
mk＝mk-1+αgk；
当|αgk|＜ε时或速度模型更新次数达到速度模型更新迭代次数N时，速度模型
更新结束；否则，进入步骤3)。
2.如权利要求1所述的一种基于最小二乘梯度更新速度模型的地震波全波形
反演方法，其特征在于：所述步骤3)中，波场误差矢量Δd的计算公式如下：
Δd＝dobs-dcal；
式中，dobs为地震观测波场数据，dcal为在迭代过程中的计算波场数据；
计算波场数据dcal的计算公式如下：
式中，为计算波场数据dcal和波场向量u之间的关系函数。
3.如权利要求1或2所述的一种基于最小二乘梯度更新速度模型的地震波全
波形反演方法，其特征在于：所述步骤5)中，地震波全波形反演的标准方程具...

【专利技术属性】
技术研发人员：杜向东，常旭，王一博，韩文明，曹向阳，张英德，张世鑫，武绍江，
申请(专利权)人：中国海洋石油总公司，中海油研究总院，
类型：发明
国别省市：北京;11