[发明专利]一种基于声波方程的线性优化隐式时空域有限差分数值模拟方法

摘要

本发明公开了一种基于声波方程的线性优化隐式时空域有限差分数值模拟方法，包括以下步骤：(1)读取参数；(2)基于菱形差分算子和二阶时间中心差分，得到时间导数的时间高阶离散格式；(3)采用空间隐式离散格式求解二阶空间导数，基于空间频散关系，优化方法求解隐式差分系数；(4)基于时间高阶离散格式和空间隐式离散格式，得到时空域频散关系；(5)采用线性优化算法对时间高阶离散格式中的差分系数进行求解；(6)利用求取的差分系数，采用混合吸收边界条件对边界反射进行吸收，按照波动方程进行递推，得到任意时刻的波场及整个地震记录；(7)记录波场快照，输出地震记录并结束。本发明的有限差分数值模拟方法精度高，稳定性好，频散误差小。

主权项

1.一种基于声波方程的线性优化隐式时空域有限差分数值模拟方法，其特征在于，包括以下几个步骤：(1)读取参数，所述参数包括正演模拟所需要的速度模型参数文件、有限差分算子长度、子波函数及子波主频、正演所采用的时间与空间步长及地震记录时长；(2)基于菱形差分算子和二阶时间中心差分，得到时间导数的时间高阶离散格式；(3)求解二阶空间导数隐式差分格式，在显式高阶差分的基础上，分母中引入二阶中心差分格式，构造空间隐式差分格式求解二阶空间导数，用于减小有限差分算子长度；进一步，基于最小二乘法，以L2范数目标函数极小化空间域频散关系，求解优化的隐式差分系数；(4)基于所述时间高阶离散格式和空间隐式离散格式，带入声波方程，构建同时具有时间高阶和空间隐式的差分递推格式，得到时空域频散关系；(5)固定空间差分系数，采用优化策略极小化整个递推格式频散关系，建立关于时间差分格式中差分系数的L2范数目标函数，利用线性优化算法对时间差分系数进行求解；(6)采用混合吸收边界条件对边界反射能量进行吸收，按照时间高阶和空间隐式差分格式离散波动方程并进行递推，得到任意时刻的波场及整个地震记录；(7)记录波场快照，输出地震记录并结束；步骤(2)中，求解二阶时间导数时采用中心差分，形式如下：其中，为压力场，τ为时间步长，t为时间；为了进一步提高时间离散精度，将菱形差分算子引入到公式(1)中，得出如下的时间高阶离散公式：其中，h为空间采样步长，v为速度，a_m,n为有限差分系数；N为菱形算子中的差分算子长度；n和m为菱形算子中的坐标位置序列；步骤(3)中，采用如下空间隐式离散格式对二阶空间导数进行近似，形式如下：其中，a′₀、a′_m和b为隐式差分系数，当b＝0时，公式(3)退化为显式差分格式；h为空间步长，x为方向坐标，M为差分算子长度；公式(3)对应的空间频散关系如下：其中，β＝k_xh，k_x为水平方向波数；极小化公式(4)对应的空间频散关系，目标函数为关于差分系数的线性凸函数，最优差分系数通过求解如下公式得到：其中，n＝1,2,L,M+1，d＝[d₁,d₂,L,d_M,d_M+1]^T＝[a′₁,a′₂,L,a′_M‑1,a′_M,b]^T，q为介于0和1之间的自然数，控制着积分范围，步骤(4)中，在求解时间导数时，采用公式(2)计算；求解空间导数时，采用公式(3)计算，将公式(2)和公式(3)带入到声波方程并进行平面波理论分析，过程如下：均匀模型中，声波方程压力P满足如下平面波方程，其中，ω是角频率，是虚数单位，k_z为垂直方向波数；m、n和l为位置和时间序列；将公式(7)带入公式(2)并进行化简得到如下公式：同样，对空间导数公式进行分析，得到如下频散关系：将上述(8)、(9)和(10)带入声波方程并简化和重组得到如下时空域频散关系，其中，ω为角频率，k_x和k_z分别为沿x和z方向的波数，r为库朗数，r＝vτ/h；步骤(5)中，首先，差分系数a′_m和b由公式(5)计算得到，在计算差分系数a_m,0和a_m,n时保持固定，此时优化差分系数a_m,0和a_m,n是线性优化问题，直接线性求解，定义如下目标函数：其中，d＝[d₁,L,d_N+(N‑1)N/2]^T＝[a_1,0,L,a_N,0,a_1,1,L,a_1,N‑1,L,a_N‑1,1]^T，表示差分系数；θ为传播角度，q为介于0和1间的实数，波数k_x和k_z与波数k之间满足如下关系：k_x＝kcos(θ),k_z＝ksin(θ),  (13)ψ(β,θ,r)及定义如下：目标函数是关于差分系数a_m,0和a_m,n的线性函数，差分系数a_m,0和a_m,n通过求解如下公式得到，