[发明专利]一种基于弥散黏滞性波动方程的FCT‑FDM正演模拟方法

权利要求书

1.一种基于弥散黏滞性波动方程的FCT-FDM正演模拟方法，其特征在于，包括以下步骤：

1)根据实际地质背景条件、岩石物理测试及其测井资料建立介质模型的参数；

2)对所设计的介质模型进行网格离散；

3)震源函数在空间上采用高斯函数，时间上采用Ricker子波，形式为：

s(x,z,t)＝g(x,z)·f(t)

式中：

f(t)＝(1-2(πf₀t)²)exp(-(πf₀t)²)

g(x,z)=exp{-[(x-x0)2+(z-z0)2]}β2]]>

式中：f₀表示Ricker子波的中心频率，模型计算中f₀＝30Hz；β为常数；(x₀,z₀)表示震源的空间位置；

4)将二维弥散黏滞性波动方程中的微分用差分近似替代，得到其相应的有限差分格式；其中空间采样步长和时间采样步长必须满足以下数值格式的稳定性条件：

0≤Δt≤min(64hυ,h2γh2+8η)]]>

其中，h为空间网格步长；

5)对二维弥散黏滞性波动方程进行FCT-FDM计算，得到每一时刻的波场值，其具体方法如下：

首先，给出弥散黏滞性波动方程的数学描述，二维弥散黏滞性波动方程其形式为：

∂2u∂t2+γ∂u∂t-η(∂3u∂x2∂t+∂3u∂z2∂t)-υ2(∂2u∂x2+∂2u∂z2)=0---(1)]]>

式中：u为波场函数；γ，η分别为弥散和黏滞性衰减系数，它们与岩石的孔隙度、渗透率以及流体的密度、粘度等有关；υ为非频散介质中波的传播速度；x，t分别为空间和时间变量；

式(1)中左端第一项表示惯性项，第二项为弥散耗损力即扩散项，第三项表示黏滞性阻尼，第四项为波动方程的弹性部分；

弥散黏滞性波动方程FCT-FDM方法的具体的计算步骤如下：

5-1)定义式(1)的数值解为给定初值为震源函数；其中，表示第n时间步在网格点(x_j,z_k)处的波场值；n为时间采样点，j为空间x方向的采样点，k为空间z方向的采样点；

5-2)利用弥散黏滞性波动方程有限差分格式计算

uj,kn+1=[(2-4a-γ(Δt)-4b]uj,kn-[1-γ(Δt)-4a]uj,kn-1-a(uj+1,kn-1+uj-1,kn-1+uj,k+1n-1+uj,k-1n-1)+(a+b)[(uj+1,kn+uj-1,kn+uj,k+1n+uj,k-1n)]---(2)]]>

式中：h为空间网格步长；Δt为时间步长；

5-3)计算漫射通量：

a)计算第n时间步的漫射通量

Pj,kn=λ1(uj+1,kn-uj,kn-uj+1,kn-1+uj,kn-1),0≤λ1≤1Qj,kn=λ1(uj,k+1n-uj,kn-uj,k+1n-1+uj,kn-1),0≤λ1≤1---(3)]]>

其中，λ₁为参数；

b)利用漫射通量平滑差分方程的解

u~j,kn+1=uj,kn+1+(Pj+12,kn-Pj-12,kn)+(Qj,k+12n-Qj,k-12n)---(4)]]>

其中，为平滑以后的解；

5-4)计算反漫射通量：

a)计算第n+1时间步的漫射通量

P~j,kn+1=λ2(uj+1,kn+1-uj,kn+1-uj+1,kn+uj,kn),0≤λ1≤1Q~j,kn+1=λ2(uj,k+1n+1-uj,kn+1-uj,k+1n+uj,kn),0≤λ2≤1---(5)]]>

式中：λ₂与λ₁类似，为参数；

b)利用和计算反漫射通量

Xj+1,kn+1=u~j+1,kn+1-uj,kn+1-(u~j+1,kn-uj,kn)Yj,k+1n+1=u~j,k+1n+1-uj,kn+1-(u~j,k+1n-uj,kn)---(6)]]>

c)利用反漫射通量修正平滑以后的解

Uj,kn+1=u~j,kn+1-(X~j+12,kn-X~j-12,kn)-(Y~j,k+12n-Y~j,k-12n)---(7)]]>

式中：为修正以后的解，且

X~j+12,kn=(sx·max{0,min[sx·Xj-12,kn,|P~j+12,kn+1|,sx·Xj+32,kn]})Y~j,k+12n=(sy·max{0,min[sy·Yj,k-12n,|Q~j,k+12n+1|,sy·Yj,k+32n]})---(8)]]>

Xj-12,kn=12(Xj-1,kn+Xj,kn),Xj+32,kn=12(Xj+1,kn+Xj+2,kn),Yj,k-12n=12(Yj,k-1n+Yj,kn),Yj,k+32n=12(Yj,k+1n+Yj,k+2n)---(9)]]>

P~j+12,kn+1=12(P~j+1,kn+1+P~j,kn+1),Q~j,k+12n+1=12(Q~j,k+1n+1+Q~j,kn+1)---(10)]]>

sx=sign(P~j+12,kn+1),sy=sign(Q~j,k+12n+1)---(11)]]>

(11)式中“sign”为符号函数；(8)式中max(·)，min(·)分别表示物理量的最大值和最小值；

5-5)重复步骤5-2)～5-4)，直到最大的计算时间步，得到整个计算区域内不同时间的波场值；

6)在计算区域的边界处使用吸收边界条件并对其进行处理，模拟实际地下介质中波的传播。

2.根据权利要求1所述的基于弥散黏滞性波动方程的FCT-FDM正演模拟方法，其特征在于：所述步骤2)中，采用矩形网格对模型进行离散。