[发明专利]一种基于弥散黏滞性波动方程的吸收边界反射方法

摘要

本发明公开了一种基于弥散黏滞性波动方程的吸收边界反射的方法，为利用该方程进行复杂介质中地震波数值模拟提供了一种处理人工边界反射波的有效途径。首先从理论上阐述了一维和二维情况下该方法用于处理人工边界反射的原理；其次，在均匀介质和复杂储层介质模型中利用有限差分法对弥散黏滞性波动方程进行数值模拟，并采用所提的吸收边界反射方法进行处理，数值模拟结果验证了该方法的有效性。该方法也适用于声波方程和Stokes方程。此外，该技术方案易于实现，可操作性强。

主权项

一种基于弥散黏滞性波动方程的吸收边界反射方法，其特征在于，包括以下步骤：1)根据求解一维弥散性黏滞性波动方程的吸收边界控制方程，推广确定二维弥散性黏滞性波动方程的吸收边界反射控制方程；其中，一维弥散性黏滞性波动方程为：∂2u∂t2+γ∂u∂t-η∂3u∂x2∂t-υ2∂2u∂x2=0---(1)]]>二维弥散性黏滞性波动方程为：∂2u∂t2+γ∂u∂t-η(∂3u∂x2∂t+∂3u∂z2∂t)-υ2(∂2u∂x2+∂2u∂z2)=0---(2)]]>式中：u为波场函数；γ，η分别为弥散和黏滞性衰减系数，它们是岩石的孔隙度、渗透率以及流体的密度、粘度等参数的函数；υ为非频散介质中波的传播速度；x，t分别为空间和时间变量；式(1)和式(2)中左端第一项表示惯性项，第二项为弥散耗损力即扩散项，第三项表示黏滞性阻尼，第四项为波动方程的弹性部分；二维弥散性黏滞性波动方程的吸收边界控制方程为：∂Sx*∂t+(d1(x)+γ)Sx*+γd1(x)sgn(t)2*Sx*-η∂Rx*∂x-υ2∂Tx*∂x=0∂Sz*∂t+(d2(z)+γ)Sz*+γd2(z)sgn(t)2*Sz*-η∂Rz*∂z-υ2∂Tz*∂z=0∂ux*∂t=Sx*,∂uz*∂t=Sz*∂ux*∂x+∂uz*∂x=Tx*+d1(x)sgn(t)2*Tx*∂ux*∂z+∂uz*∂z=Tz*+d2(z)sgn(t)2*Tz*∂Sx*∂x+∂Sz*∂x=Rs*+d1(x)sgn(t)2*Rx*∂Sx*∂z+∂Sz*∂z=Rz*+d2(z)sgn(t)2*Tz*---(3)]]>式中：d1(x)，d2(z)分别为沿x，z方向的衰减函数；2)建立介质模型，在其周围加入吸收边界层，吸收边界层的厚度取值至少为一个波长；3)选择吸收边界中的吸收衰减函数，其中吸收衰减函数包括指数函数、对数函数、正/余弦函数或其中的一种或几种的复合函数；4)对所计算的介质模型进行网格离散；5)利用有限差分法对二维弥散黏滞性波动方程进行数值模拟，在计算区域内部利用式(4)弥散黏滞性波动方程的有限差分格式计算，在吸收边界层内利用式(3)计算，以消除人工边界反射波；其中，式(4)如下：uj,mn+1=[2-4a-γ(Δt)-4b]uj,mn-[1-γ(Δt)-4a]uj,mn-1-a(uj+1,mn-1+uj-1,mn-1+uj,m+1n-1+uj,m-1n-1)+(a+b)[(uj+1,mn+uj-1,mn+uj,m+1n+uj,m-1n)]---(4)]]>式中：表示第n时间步在网格点(xj,zm)处的波场值；n为时间采样点，j为空间x方向的采样点，m为空间z方向的采样点；震源函数采用爆炸源，在空间上采用高斯函数，时间上采用Ricker子波，震源函数的形式为：s(x,z,t)＝g(x,z)·f(t)式中：f(t)＝(1?2(πf0t)2)exp(?(πf0t)2)g(x,z)=exp[-[(x-x0)2+(z-z0)2]}β2]]>式中：f0表示Ricker子波的中心频率，模型计算中β为常数；(x0,z0)表示震源的空间位置。