【技术实现步骤摘要】
基于有限元的逆时偏移方法和装置
本专利技术属于油气勘探开发中的地震成像领域，更具体地，涉及一种基于有限元的逆时偏移方法和一种基于有限元的逆时偏移装置。
技术介绍
目前已经工业化应用的地震叠前逆时偏移方法均是基于不同差分格式的有限差分法。有限差分逆时偏移成像方法计算简单，效率较高，满足工业生产时的时效性需求。但是，这种有限差分逆时偏移成像方法在面对极不规则地表或地下岩体时，成像精度不高。也有少数研究人员在研究基于有限元的叠前逆时偏移方法，其边界适应性强，可以对极不规则地表或地下岩体进行精确成像。但现有的基于有限元的叠前逆时偏移方法计算量巨大，并行效率低，无法满足工业应用对时效性要求，还停留在针对二维模型的实验研究阶段。
技术实现思路
有鉴于此，本申请提出了一种基于有限元的可在工业应用中推广的高效逆时偏移成像方法，本申请还提出了相应的装置。根据本申请的一个方面，提出了一种基于有限元的逆时偏移方法，所述方法包括：将连续的三维空间分割成相互连通但不重合的有限个单元，并建立空间离散的有限元控制方程；根据时间积分法外推有限元控制方程，得到空间和时间均离散的线性方程组；通过分割面将所述三维空间划分为多个分块，并将每个分块中各个单元节点的位移更新值用q个位移增量来表示；根据最小势能原理对所述线性方程组进行变换，并代入用q个位移增量来表示单元节点的位移更新值，得到多级并行线性方程组：其中，Set(I)是第I分块中所有单元节点的集合，Set(J)是第J分
【技术保护点】
1.一种基于有限元的逆时偏移方法，其特征在于，所述方法包括：/n将连续的三维空间分割成相互连通但不重合的有限个单元，并建立空间离散的有限元控制方程；/n根据时间积分法外推有限元控制方程，得到空间和时间均离散的线性方程组；/n通过分割面将所述三维空间划分为多个分块，并将每个分块中各个单元节点的位移更新值用q个位移增量来表示；/n根据最小势能原理对所述线性方程组进行变换，并代入用q个位移增量来表示单元节点的位移更新值，得到多级并行线性方程组：/n
【技术特征摘要】
1.一种基于有限元的逆时偏移方法，其特征在于，所述方法包括：
将连续的三维空间分割成相互连通但不重合的有限个单元，并建立空间离散的有限元控制方程；
根据时间积分法外推有限元控制方程，得到空间和时间均离散的线性方程组；
通过分割面将所述三维空间划分为多个分块，并将每个分块中各个单元节点的位移更新值用q个位移增量来表示；
根据最小势能原理对所述线性方程组进行变换，并代入用q个位移增量来表示单元节点的位移更新值，得到多级并行线性方程组：
其中，Set(I)是第I分块中所有单元节点的集合，Set(J)是第J分块中所有单元节点的集合，为第I分块中单元节点i的第m个位移增量模式，为第J分块中单元节点j的第l个位移增量模式，kij为对应的整体刚度，为第J分块的第l个位移增量的系数，si为震源在单元节点i上施加的压力，为位移更新前单元节点i的位移近似值，B为分块的总数；
调用多个从处理器分别根据下式处理各个分块的数据，每个从处理器对应一个分块：
调用主控处理器接收各个从处理器的处理结果，以求解多级并行线性方程组。
2.根据权利要求1所述的方法，其特征在于，所述分割面被置于单元内。
3.根据权利要求1所述的方法，其特征在于，所述将连续的三维空间分割成相互连通但不重合的有限个单元，并建立空间离散的有限元控制方程，包括：
得到时空域中直角坐标系下的声波方程：
其中，P是声压，v是波速，t是时间变量，x、y、z是空间变量，s是震源函数；
采用部分离散，构造如下近似解：
其中，Nj(x,y,z)是单元节点j上的插值函数，以下简写为Nj，dj(t)是单元节点j的位移，以下简写为dj，nd是单元节点总数；
将公式2代入公式1，产生余量R：
根据加权余量法，利用近似解的插值函数作权函数，令余量R在求解区域Ω上的加权积分为零，即
∫ΩNiRdΩ＝0(i＝1,2,…,nd)公式4
将公式3式代入公式4，得：
对公式5等号左端中间三项作分部积分，得
其中，nx、ny、nz为边界外法线的方向余弦，Γ是Ω的外边界；
将公式2和公式6式代入公式5，得到空间离散的有限元控制方程：
其中，为dj的二阶偏导。
4.根据权利要求3所述的方法，其特征在于，所述根据时间积分法外推有限元控制方程，得到空间和时间离散的线性方程组，包括：
将有限元控制方程改写为矩阵形式：
其中，
其中，N是各个单元节点的插值函数矩阵，d是各个单元节点的位移矩阵，是各个单元节点的位移二阶偏导矩阵，ρ是密度，是哈密顿算子，Ωe表示积分区域是一个单元；
在沿时间轴的计算中，使用隐式NewMark方法实现时间积分，在时间区域内NewMark方法使用以下假设：
根据公式10进一步可以得到
其中，是各个单元节点的位移一阶偏导矩阵，Δt表示时间增量，下标t和t+Δt表示时间变量的取值，α和δ为常数系数；
根据公式10和公式11整理公式8和公式9，得到：
其中有：
求解公式12后更新得到：
其中，为整体刚度矩阵，c0、c1、c2、c3、c4、c5、c6、c7均为常数系数，C为阻尼矩阵；
由于实际使用随机速度边界条件，因而计算中不出现阻尼矩阵C，公式12简化为：
边界条件为
其中，d0表示各个节点的位移的初始值矩阵，表示各个节点位移二阶偏导的初始值矩阵。
5.根据权利要求1所述的方法，其特征在于，所述将每个分块中各个单元节点的位移更新值用q个位移增量来表示，包括：
用下式表示第I分块中单元内节点i的位移di：
其中，为第I分块中第m个位移增量的系数，为第I分块中节点i的第m个位移增量模式。
6....
【专利技术属性】
技术研发人员：刘立民，刘定进，姚晓龙，李博，
申请(专利权)人：中国石油化工股份有限公司，中国石油化工股份有限公司石油物探技术研究院，
类型：发明
国别省市：北京;11
还没有人留言评论。发表了对其他浏览者有用的留言会获得科技券。