基于3D Gauss-FFT的三维重力场正演方法技术

本发明专利技术提供一种基于3D Gauss‑FFT的三维重力场正演方法，包括以下步骤：步骤S100：三维场源体建模；步骤S200:3D Gauss‑FFT正变换计算三维重力场偏移频谱；步骤S300：3D Gauss‑FFT逆变换计算单个高斯节点的三维重力场；步骤S400:累加求和得到三维重力场分布。该方法首先提出3D Gauss‑FFT算法，即在每一个积分区域用高斯型求积公式代替传统的矩形求积公式，大大提高三维重力场计算精度，同时由于FFT算法内核的使用，计算效率得到保证，模型试验表明本方法正演精度与空间域解析解相当。

【技术实现步骤摘要】
基于3DGauss-FFT的三维重力场正演方法
本专利技术涉及勘探地球物理
，具体的涉及一种基于3DGauss-FFT的三维重力场正演方法。
技术介绍
重力勘探作为一项最基本的地球物理方法，已被广泛运用于固体矿产资源勘探、岩石圈壳幔结构研究、地形校正、大地水准面测量、水文地质勘察、考古勘察等。测量得到的重力场数据与地下密度异常体有着密切关系，可被直接用于反演得到地下三维密度分布，从而推断地下构造运动、岩石圈结构构造等科学问题。正演的速度和精度决定了反演的可行性及可信度，因此在一个高效的反演程序中，快速高精度正演算法尤为重要。在直角坐标系下，地下场源场被剖分为规则的直立棱柱体，传统的重力场正演方法有空间域解析解法，频率域FFT法。然而空间域解析解法计算零误差，但计算异常耗时；频率域FFT法计算效率高而计算精度难以保证。这一矛盾在三维重力场反演中尤为显著，不精确的正演结果以及低效率的正演方法将制约着反演的进行。在当前计算机硬件难以突破的背景下，想要大幅度提高计算效率，只能寻求新的正演方法。
技术实现思路
本专利技术的目的在于提供一种基于Gauss-FFT的三维重力场正演方法，该专利技术解决了现有直角坐标系下三维重力场正演算法正演效率低，计算精度差的技术问题。本专利技术提供一种基于Gauss-FFT的三维重力场正演方法，包括以下步骤：步骤S100：根据地下三维异常体的形态和尺寸，设置正演模型并剖分所述正演模型得到等体积的多个直立棱柱体小单元体，以任一所述小单元体的几何中心处为高斯节点，设定所述正演模型中高斯节点个数，每一个方向上高斯节点数为2，可以更多，查询[-1，1]上高斯节点系数表，将其转换为[0，1]上相应的节点值和系数值；步骤S200：对第i个高斯节点采用3DGauss-FFT正变换计算三维重力场偏移频谱：步骤S210：分别设置X、Y、Z方向上第i个高斯节点和第i个高斯节点处的高斯系数数值对(λix,αix),(λiy,αiy),(λiz,αiz)，其中，λix为X方向上第i个高斯节点处的高斯系数数，λiy为Y方向上第i个高斯节点处的高斯系数数，λiz为Z方向上第i个高斯节点处的高斯系数数，αix为X方向上第i个高斯节点的值，αiy为Y方向上第i个高斯节点的值，αiz为Z方向上第i个高斯节点的值，离散密度ρ(xm,yn,zl)乘以高斯偏移因子得到第i个高斯点上的偏移密度其中，e为自然底数；xm为X方向上第m段的中心点坐标，yn为Y方向上第n段的中心点坐标，zl为Z方向上第l段的中心点坐标，αix为X方向上第i个高斯节点的值，αiy为Y方向上第i个高斯节点的值，αiz为Z方向上第i个高斯节点的值；Δkx为kx方向上频率的间隔，Δky为ky方向上的频率间隔，Δkz为kz方向上频率的间隔。步骤S220：对所述偏移密度进行三维离散傅里叶变换(3DDFT)，得到第i个高斯点上的密度分布的偏移频谱步骤S230：对于所述密度分布的偏移频谱乘以大地滤波因子，得到第i个高斯点上的三维重力异常偏移频谱其中，kxp为kx方向上第p个点处的频率值，kyq为ky方向上第q个点处的频率值，kzw为kz方向上第w个点处的频率值；步骤S300：对第i个高斯节点采用3DGauss-FFT逆变换计算第i个高斯节点上的重力场响应；步骤S400：取i＝i+1重复所述步骤S200～S300，直至i＝n时停止，其中n为总的高斯节点数，累加所得各高斯节点上的重力场响应得到所述地下三维异常体的三维重力场分布。进一步地，所述步骤S100包括以下步骤：步骤S110：根据地下三维异常体的形态和尺寸，设置用于包含所述三维异常体的棱柱体目标区域，所设置的棱柱体目标区域能将整个三维异常体完全嵌入其中；步骤S120：设置所述三维异常体的剖分段数，确定所述三维异常体X、Y、Z方向上的剖分间隔，根据所述剖分段数和所述剖分间隔将所述目标区域剖分成体积相等的直立棱柱体小单元体，每一个小单元体的几何中心处为高斯节点，剖分段数的设置可根据问题的实际需要以及计算机的实际性能进行设置；步骤S130：根据所述三维异常体的密度分布，对每一个剖分后的小直立棱柱体单元密度进行赋值，得到ρ(xm,yn,zl)。进一步地，所述步骤S300包括以下步骤：步骤S310：对第i个高斯节点上的三维重力场偏移频谱做三维离散傅里叶逆变换(3DIDFT)，得到第i个高斯点处三维重力场偏移频谱空间域值；步骤S320：将所述第i个高斯点处三维重力场偏移频谱空间域值乘以高斯逆偏移因子和高斯系数，得到第i个高斯节点上的重力场响应。进一步地，所述步骤S220中所述三维离散傅里叶变换(3DDFT)包括以下步骤：步骤S221：将排列好的偏移密度数据输入到Fortran库函数fftw3，程序将自动得到偏移密度的频谱。进一步地，所述步骤S300中3DIDFT逆变换包括以下步骤：步骤S310：将排列好的第i个点处所对应的三维重力异常频谱输入到Fortran库函数fftw3中，程序将自动计算其逆变换值并输出；步骤S320：将程序输出的结果再乘以第i个高斯点所对应的高斯系数以及逆偏移因子，得到第i个高斯节点上的重力异常分布。相对现有技术，本专利技术的技术效果：本专利技术提供的基于Gauss-FFT的三维重力场正演方法，首先提出3DGauss-FFT算法，即在每一个积分区域用高斯型求积公式代替传统的矩形求积公式，大大提高三维重力场计算精度，同时由于FFT算法内核的使用，计算效率得到保证，模型试验表明本方法正演精度与空间域解析解相当。本专利技术提供的基于Gauss-FFT的三维重力场正演方法，核心内容为：针对传统基于FFT的三维重力场正演方法，提出用3DGauss-FFT算法代替FFT算法，在每一个单元体上用更高精度的高斯型积分代替传统FFT算法中的矩形积分，有效减小了传统离散傅里叶变换中由于边界效应、强加周期、截断效应等所引起的诡源效应，大幅度提高了正演计算精度，同时由于FFT内核的使用，正演计算效率也得到了保证。模型试验表明，本专利技术所提出的方法使得三维重力场正演精度相比于传统FFT算法提高了近两个数量级，有效压制了诡源效应。本专利技术所提出的方法可直接用于三维重力场正演，井地重力场联合反演。本专利技术提供的基于Gauss-FFT的三维重力场正演方法，在传统重力场频域正演中，使用快速傅里叶变换来提高计算效率。在现有技术中采用传统离散傅里叶变换算法时，在每一个积分区域中使用矩形求积公式，因此其正演精度较低，尤其是在截断边界处。本专利技术提出了3DGauss-FFT方法，该方法在每一个积分区域中使用高斯型求积公式代替传统的矩形求积公式，使得计算精度大大提高。本专利技术提供的基于Gauss-FFT的三维重力场正演方法，由于传统FFT内核的调用，本专利技术所提出的方法同样具有传统频域正演算法的高效性。模型算例表明，本专利技术所提出的方法能够有效减小诡源效应的影响，使得重力异常分布能够正确反映地下场源分布。诡源效应的消除对于三维重力异常反演至关重要，正确的异常分布才有可能得到正确的密度异常体分布。此外，与传统重力场频域算法相比，本专利技术所提出的方法在正演精度方面提高了两个数量级，在计算效率方法与传统方法相当，算例证明本专利技术所提出的方法是一种快速高精度的三维重力场正演方法，可本文档来自技高网...

【技术保护点】
1.一种基于3D Gauss‑FFT的三维重力场正演方法，其特征在于，包括以下步骤：步骤S100：根据地下三维异常体的形态和尺寸，设置正演模型并剖分所述正演模型得到等体积的多个直立棱柱体小单元体，以任一所述小单元体的几何中心处为高斯节点，设定所述正演模型中高斯节点个数，每一个方向上高斯节点数为2，可以更多，查询[‑1，1]上高斯节点系数表，将其转换为[0，1]上相应的节点值和系数值；步骤S200：对第i个高斯节点采用3D Gauss‑FFT正变换计算三维重力场偏移频谱：步骤S210：分别设置X、Y、Z方向上第i个高斯节点和第i个高斯节点处的高斯系数数值对(λix,αix),(λiy,αiy),(λiz,αiz)，其中，λix为X方向上第i个高斯节点处的高斯系数数，λiy为Y方向上第i个高斯节点处的高斯系数数，λiz为Z方向上第i个高斯节点处的高斯系数数，αix为X方向上第i个高斯节点的值，αiy为Y方向上第i个高斯节点的值，αiz为Z方向上第i个高斯节点的值，离散密度ρ(xm,yn,zl)乘以高斯偏移因子

【技术特征摘要】
1.一种基于3DGauss-FFT的三维重力场正演方法，其特征在于，包括以下步骤：步骤S100：根据地下三维异常体的形态和尺寸，设置正演模型并剖分所述正演模型得到等体积的多个直立棱柱体小单元体，以任一所述小单元体的几何中心处为高斯节点，设定所述正演模型中高斯节点个数，每一个方向上高斯节点数为2，可以更多，查询[-1，1]上高斯节点系数表，将其转换为[0，1]上相应的节点值和系数值；步骤S200：对第i个高斯节点采用3DGauss-FFT正变换计算三维重力场偏移频谱：步骤S210：分别设置X、Y、Z方向上第i个高斯节点和第i个高斯节点处的高斯系数数值对(λix,αix),(λiy,αiy),(λiz,αiz)，其中，λix为X方向上第i个高斯节点处的高斯系数数，λiy为Y方向上第i个高斯节点处的高斯系数数，λiz为Z方向上第i个高斯节点处的高斯系数数，αix为X方向上第i个高斯节点的值，αiy为Y方向上第i个高斯节点的值，αiz为Z方向上第i个高斯节点的值，离散密度ρ(xm,yn,zl)乘以高斯偏移因子得到第i个高斯点上的偏移密度其中j为虚数单位；步骤S220：对所述偏移密度进行三维离散傅里叶变换(3DDFT)，得到第i个高斯点上的密度分布的偏移频谱步骤S230：对于所述密度分布的偏移频谱乘以大地滤波因子，得到第i个高斯点上的三维重力异常偏移频谱步骤S300：对第i个高斯节点采用3DGauss-FFT逆变换计算第i个高斯节点上的重力场响应；步骤S400：取i＝i+1重复所述步骤S200～S300，直至i＝n时停止，其中n为总的高斯节点数，累加所得各高斯节点上的重力场响应得到所述地下三维异常体的三维重力场分布。2.根据权利要求1所述的基于3DGauss-FFT的三维重力场正演方法，其特征在于，所述步骤S100包括以下步骤：步骤S110：根据地下三维异常体的形态和尺寸...

【专利技术属性】
技术研发人员：赵广东，柳建新，
申请(专利权)人：中南大学，
类型：发明
国别省市：湖南,43