正交介质地震波场模拟频散压制方法与流程

技术特征：

1.正交介质地震波场模拟频散压制方法，其特征在于，该正交介质地震波场模拟频散压制方法包括：

步骤1，利用岩石物理理论和裂隙等效理论将裂隙物性参数转化为正交介质参数，为正交介质正演模拟提供模型；

步骤2，根据地震波动力学理论结合裂隙介质刚度矩阵推导正交介质弹性波一阶速度-应力方程；

步骤3，利用有限差分法求解正交介质弹性波方程模拟地震波在裂隙介质中的传播过程；

步骤4，利用最小二乘优化方法优化正交介质交错网格差分系数从而压制频散提高模拟精度；

步骤5，进行模拟测试。

2.根据权利要求1所述的正交介质地震波场模拟频散压制方法，其特征在于，在步骤1中，借助线性滑动理论给出了各向同性介质背景下的两组正交直立裂隙的弹性系数矩阵：

$C=\left[\begin{array}{cccccc}{C}_{11}& C_{12}& C_{13}& & & \\ C_{12}& C_{22}& C_{23}& & & \\ C_{13}& C_{23}& C_{33}& & & \\ & & & C_{44}& & \\ & & & & C_{55}& \\ & & & & & C_{66}\end{array}\right]=\left(\begin{array}{cc}{\tilde{C}}_1& 0\\ 0& {\tilde{C}}_2\end{array}\right)---\left(1\right)$

其中C为弹性系数矩阵，C₁₁,C₂₂,C₃₃,C₄₄,C₅₅,C₆₆,C₁₂,C₁₃,C₂₃为弹性系数矩阵元素，为子矩阵，其又可表示为：

${\tilde{C}}_1=\frac{1}{d}\left[\begin{array}{ccc}(\mathrm{\λ}+2\mathrm{\μ})l_1{m}_3& {\mathrm{\λl}}_1{m}_1& {\mathrm{\λl}}_1{m}_2\\ {\mathrm{\λl}}_1{m}_1& (\mathrm{\λ}+2\mathrm{\μ})l_3{m}_1& {\mathrm{\λl}}_2{m}_1\\ {\mathrm{\λl}}_1{m}_2& {\mathrm{\λl}}_2{m}_1& (\mathrm{\λ}+2\mathrm{\μ})(l_1{m}_3-l_4)\end{array}\right]---\left(2\right)$

${\tilde{C}}_2=\left[\begin{array}{ccc}\mathrm{\μ}(1-{\mathrm{\Δ}}_{T2})& 0& 0\\ 0& \mathrm{\μ}(1-{\mathrm{\Δ}}_{T1})& 0\\ 0& 0& \frac{\mathrm{\μ}(1-{\mathrm{\Δ}}_{T1})(1-{\mathrm{\Δ}}_{T2})}{(1-{\mathrm{\Δ}}_{T1}{\mathrm{\Δ}}_{T2})}\end{array}\right]---\left(3\right)$

其中λ和μ分别为裂隙所在背景介质的拉梅参数，Δ_N和Δ_T分别是法向柔度和切向柔度，它们与裂隙充填物有关，变量的下标N1,N2,T1,T2分别表示第一组和第二组裂隙法向和切向，d＝1-r²Δ_N1Δ_N2,r＝1-2g,g＝μ/(λ+2μ)，l₁,l₂,l₃,l₄,m₁,m₂,m₃,d为过渡参数，其可以表示为：

$\begin{array}{c}{l}_1=1-{\mathrm{\Δ}}_{N1},l_2=1-r-{\mathrm{\Δ}}_{N1},l_3=1-r^2{\mathrm{\Δ}}_{N1},l_4=4r^2{g}^2{\mathrm{\Δ}}_{N1}{\mathrm{\Δ}}_{N2}\\ m_1=1-{\mathrm{\Δ}}_{N2},m_2=1-{\mathrm{r\Δ}}_{N2},m_3=1-r^2{\mathrm{\Δ}}_{N2}\end{array}---\left(4\right)$

切向柔量和法向柔量的表达式，其中K为体积模量：

${\mathrm{\Δ}}_{Ni}=\frac{(\mathrm{\λ}+2\mathrm{\μ})K_{Ni}}{1+(\mathrm{\λ}+2\mathrm{\μ})K_{Ni}},{\mathrm{\Δ}}_{Ti}=\frac{{\mathrm{\μK}}_{Ti}}{1+{\mathrm{\μK}}_{Ti}}---\left(5\right)$

当裂隙满足Hudson理论的假设时，可以从Hudson理论中给出Δ_T,Δ_N的表达式：

(1)当裂缝中填充较小体积模量和剪切模量的固体时：

${\mathrm{\Δ}}_N=\frac{4}{3g\[1-g+(k^{\′}+\frac{4{\mathrm{\μ}}^{\′}}{3})/\left(\mathrm{\π}d\mathrm{\μ}\right)\]}e,{\mathrm{\Δ}}_T=\frac{16}{3\[3-2g+\frac{4{\mathrm{\μ}}^{\′}}{\mathrm{\π}d\mathrm{\μ}}\]}e---\left(6\right)$

(2)当裂缝为干裂缝时：

${\mathrm{\Δ}}_N=\frac{4}{3g(1-g)}e,{\mathrm{\Δ}}_T=\frac{16}{3(3-2g)}e---\left(7\right)$

(3)当裂缝中填充无粘滞流体时：

Δ_N＝0，

其中，k'和μ'分别为填充介质的体积模量和剪切模量，e和d分别为裂缝密度和裂缝的纵横比。

3.根据权利要求1所述的正交介质地震波场模拟频散压制方法，其特征在于，在步骤2中，根据应力应变关系：

$\left[\begin{array}{c}{\mathrm{\τ}}_{xx}\\ {\mathrm{\τ}}_{yy}\\ {\mathrm{\τ}}_{zz}\\ {\mathrm{\τ}}_{yz}\\ {\mathrm{\τ}}_{xz}\\ {\mathrm{\τ}}_{xy}\end{array}\right]=\left[\begin{array}{cccccc}{c}_{11}& c_{12}& c_{13}& 0& 0& c_{16}\\ c_{12}& c_{22}& c_{23}& 0& 0& c_{26}\\ c_{13}& c_{23}& c_{33}& 0& 0& c_{36}\\ 0& 0& 0& c_{44}& c_{45}& 0\\ 0& 0& 0& c_{45}& c_{55}& 0\\ c_{16}& c_{26}& c_{36}& 0& 0& c_{66}\end{array}\right]\left[\begin{array}{c}{e}_{xx}\\ e_{yy}\\ e_{zz}\\ e_{yz}\\ e_{xz}\\ e_{xy}\end{array}\right]---\left(9\right)$

其中τ_xx,τ_yy,τ_zz,τ_xy,τ_xz,τ_yz为应力张量，e_xx,e_yy,e_zz,e_xy,e_xz,e_yz为应变张量，C₁₁,C₂₂,C₃₃,C₄₄,C₅₅,C₆₆,C₁₂,C₁₃,C₂₃为弹性系数矩阵元素，

应力位移关系：

$\mathrm{\ρ}\frac{{\∂}^2}{\∂t^2}\left[\begin{array}{c}{u}_x\\ u_y\\ u_z\end{array}\right]=\left[\begin{array}{cccccc}\frac{\∂}{\∂x}& 0& 0& 0& \frac{\∂}{\∂z}& \frac{\∂}{\∂y}\\ 0& \frac{\∂}{\∂y}& 0& \frac{\∂}{\∂z}& 0& \frac{\∂}{\∂x}\\ 0& 0& \frac{\∂}{\∂z}& \frac{\∂}{\∂y}& \frac{\∂}{\∂x}& 0\end{array}\right]\left[\begin{array}{c}{\mathrm{\τ}}_{xx}\\ {\mathrm{\τ}}_{yy}\\ {\mathrm{\τ}}_{zz}\\ {\mathrm{\τ}}_{yz}\\ {\mathrm{\τ}}_{xz}\\ {\mathrm{\τ}}_{xy}\end{array}\right]+\mathrm{\ρ}\left[\begin{array}{c}{F}_x\\ F_y\\ F_z\end{array}\right]---\left(10\right)$

其中ρ为密度，u_x,u_y,u_z为位移变量，F_x,F_y,F_z为震源变量，

位移应变关系：

$\left[\begin{array}{c}{e}_{xx}\\ e_{yy}\\ e_{zz}\\ e_{yz}\\ e_{xz}\\ e_{xy}\end{array}\right]=\left[\begin{array}{ccc}\frac{\∂}{\∂x}& 0& 0\\ 0& \frac{\∂}{\∂y}& 0\\ 0& 0& \frac{\∂}{\∂z}\\ 0& \frac{\∂}{\∂z}& \frac{\∂}{\∂y}\\ \frac{\∂}{\∂z}& 0& \frac{\∂}{\∂x}\\ \frac{\∂}{\∂y}& \frac{\∂}{\∂x}& 0\end{array}\right]\left[\begin{array}{c}{u}_x\\ u_y\\ u_z\end{array}\right]---\left(11\right)$

即可推导三维正交介质弹性波一阶速度-应力方程：

$\left\{\begin{array}{c}\mathrm{\ρ}\frac{\∂v_x}{\∂t}=\frac{\∂{\mathrm{\τ}}_{xx}}{\∂x}+\frac{\∂{\mathrm{\τ}}_{xz}}{\∂z}+\frac{\∂{\mathrm{\τ}}_{xy}}{\∂y}+{\mathrm{\ρF}}_x\\ \mathrm{\ρ}\frac{\∂v_y}{\∂t}=\frac{\∂{\mathrm{\τ}}_{yy}}{\∂y}+\frac{\∂{\mathrm{\τ}}_{yz}}{\∂z}+\frac{\∂{\mathrm{\τ}}_{xy}}{\∂x}+{\mathrm{\ρF}}_y\\ \mathrm{\ρ}\frac{\∂v_z}{\∂t}=\frac{\∂{\mathrm{\τ}}_{zz}}{\∂z}+\frac{\∂{\mathrm{\τ}}_{yz}}{\∂y}+\frac{\∂{\mathrm{\τ}}_{xz}}{\∂x}+{\mathrm{\ρF}}_z\end{array}\right.$

$\left\{\begin{array}{c}\frac{\∂{\mathrm{\τ}}_{xx}}{\∂t}=c_{11}\frac{\∂v_x}{\∂x}+c_{12}\frac{\∂v_y}{\∂y}+c_{13}\frac{\∂v_z}{\∂z}+c_{16}(\frac{\∂v_x}{\∂y}+\frac{\∂v_y}{\∂x})\\ \frac{\∂{\mathrm{\τ}}_{yy}}{\∂t}=c_{12}\frac{\∂v_x}{\∂x}+c_{22}\frac{\∂v_y}{\∂y}+c_{23}\frac{\∂v_z}{\∂z}+c_{26}(\frac{\∂v_x}{\∂y}+\frac{\∂v_y}{\∂x})\\ \frac{\∂{\mathrm{\τ}}_{zz}}{\∂t}=c_{13}\frac{\∂v_x}{\∂x}+c_{23}\frac{\∂v_y}{\∂y}+c_{33}\frac{\∂v_z}{\∂z}+c_{36}(\frac{\∂v_x}{\∂y}+\frac{\∂v_y}{\∂x})\\ \frac{\∂{\mathrm{\τ}}_{yz}}{\∂t}=c_{44}\left(\frac{\∂v_y}{\∂z}+\frac{\∂v_z}{\∂y}\right)+c_{45}(\frac{\∂v_x}{\∂z}+\frac{\∂v_z}{\∂x})\\ \frac{\∂{\mathrm{\τ}}_{xz}}{\∂t}=c_{45}\left(\frac{\∂v_y}{\∂z}+\frac{\∂v_z}{\∂y}\right)+c_{55}(\frac{\∂v_x}{\∂z}+\frac{\∂v_z}{\∂x})\\ \frac{\∂{\mathrm{\τ}}_{xy}}{\∂t}=c_{16}\frac{\∂v_x}{\∂x}+c_{26}\frac{\∂v_y}{\∂y}+c_{36}\frac{\∂v_z}{\∂z}+c_{66}(\frac{\∂v_x}{\∂y}+\frac{\∂v_y}{\∂x})\end{array}\right.$

其中，v_x,v_y,v_z为速度变量，t表示时间变量。

4.根据权利要求1所述的正交介质地震波场模拟频散压制方法，其特征在于，在步骤3中，利用泰勒展开得到时间2阶和空间2N阶差分形式：

时间2阶差分：

$\frac{\∂u\left(t\right)}{\∂t}=\frac{u(t+\mathrm{\Δ}t)-u(t-\mathrm{\Δ}t)}{2\mathrm{\Δ}t}---\left(12\right)$

空间2N阶差分：

$\frac{\∂u}{\∂x}=\frac{1}{\mathrm{\Δ}x}\underset{n=1}{\overset{N}{\Σ}}{c}_n\[u(x+(2n-1\left)\frac{\mathrm{\Δ}x}{2}\right)-u(x-(2n-1\left)\frac{\mathrm{\Δ}x}{2}\right)\]+O\left({\mathrm{\Δx}}^{2N}\right)---\left(13\right)$

其中O(Δx^2N)为泰勒展开误差项，N阶差分系数c_m可以通过以下方程求取：

$\left[\begin{array}{ccccc}{1}^1& 3^1& 5^1& \mathrm{...}& {(2N-1)}^1\\ 1^3& 3^3& 5^3& \mathrm{...}& {(2N-1)}^3\\ 1^5& 3^5& 5^5& \mathrm{...}& {(2N-1)}^5\\ \mathrm{...}& \mathrm{...}& \mathrm{...}& \mathrm{...}& \mathrm{...}\\ 1^{(2N-1)}& 3^{(2N-1)}& 5^{(2N-1)}& \mathrm{...}& {(2N-1)}^{(2N-1)}\end{array}\right]\left[\begin{array}{c}{c}_1\\ c_2\\ c_3\\ .\\ .\\ .\\ c_N\end{array}\right]=\left[\begin{array}{c}1\\ 0\\ 0\\ .\\ .\\ .\\ 0\end{array}\right]---\left(14\right).$

5.根据权利要求1所述的正交介质地震波场模拟频散压制方法，其特征在于，在步骤4中，将最小二乘优化差分系数的方法应用于三维正交介质模拟中，一阶空间导数Taylor展开：

$\frac{\∂u}{\∂x}\≈\frac{1}{h}\underset{m=1}{\overset{M}{\Σ}}{c}_m\[u(x+mh-0.5h)-u(x-mh+0.5h)\]---\left(15\right)$

其中h为空间采样间隔，c_m为空间差分系数，M为差分阶数，

将平面波解u(x)＝u₀e^ikx带入，其中_k为波数，结合欧拉公式

e^±ix＝cos x±i sin x，化简可得：

$\frac{kh}{2}\≈\underset{m=1}{\overset{M}{\Σ}}{a}_m\sin \[kh(2m-1)/2\]---\left(16\right)$

定义过渡函数β和为：

可得：

定义误差函数：

其中E表示误差，a_M为空间差分系数，根据最小二乘基本原理求解：

$\frac{\∂E}{\∂c_m}=0---\left(20\right)$

写成矩阵形式即：

其中c_m为空间差分系数，M为差分阶数，算法和分别为：