矩形棱柱体重力梯度张量异常正演计算公式

引用本文

舒晴, 朱晓颖, 周坚鑫, 高维, 尹航. 矩形棱柱体重力梯度张量异常正演计算公式[J]. 物探与化探, 2017,39(6): 1217-1222.
SHU Qing, ZHU Xiao-Ying, ZHOU Jian-Xin, GAO Wei, YIN Hang. Forward modeling expressions for right rectangular prism with constant density[J]. Geophysical and Geochemical Exploration, 2017,39(6): 1217-1222. 复制到剪切板

Doi:10.11720/wtyht.2015.6.20
Permissions

《物探与化探》编辑部

矩形棱柱体重力梯度张量异常正演计算公式

舒晴^1,², 朱晓颖², 周坚鑫², 高维², 尹航²

1.吉林大学地球探测科学与技术学院,吉林长春 130021

2.中国国土资源航空物探遥感中心,北京 100083

作者简介: 舒晴(1979-),男,高级工程师,吉林大学在读博士,2004年中国地质大学(北京)地球物理与信息技术学院硕士毕业,长期从事航空地球物理方法技术研究与地质解释工作。

收稿日期: 2015-06-08

基金: 国家高技术研究发展计划(2011AA060501)、中国地质调查局地质矿产调查评价项目(1212011087013)

摘要

矩形棱柱体重力梯度张量异常常用的正演计算公式含有解析奇点。结合理论推导过程详细分析了“奇点”存在的原因,并基于前人提出的长方体Δ T场及其梯度场无解析奇点理论表达式,建立了矩形棱柱体重力梯度张量无解析“奇点”公式,通过模型计算结果对比和理论分析证明了其正确性。

关键词: 重力梯度张量异常; 矩形棱柱体; 解析奇点

中图分类号:P632 文献标志码:A 文章编号:1000-8918(2015)06-1217-06

Forward modeling expressions for right rectangular prism with constant density

SHU Qing^1,², ZHU Xiao-Ying², ZHOU Jian-Xin², GAO Wei², YIN Hang²

1.College of Geoexploration Science and Technology,Jilin University,Changchun 130021,China

2.China Aero Geophysics Survey and Remote Sensing Center for Land and Resources,Beijing 100083,China

Abstract

General forward modeling expressions for right rectangular prism are derived from Newton's law of universal gravitation.In this paper, the authors dealt with reasons which lead to analytic singular points in the existing formulae and established non-singularity expressions based on the special integration method presented by Guo Zhihong.The validity of the new expressions was tested by model comparative study and theoretical analysis.

Keyword: gravity gradient tensor anomaly; rectangular prism; analytic singular point

Show Figures

矩形棱柱体作为最常用的模型构建单元被广泛应用于重力梯度正演数值模拟之中。笔者研究发现, 目前被广泛引用的矩形棱柱体重力梯度张量异常计算公式^{[1, 2, 3, 4, 5]}均源自Donald Plouff^[6], 该公式虽然简洁、对称, 但文献中仅给出了结果, 没有详细地推导过程, 且存在解析奇点^[7], 在实际应用中需要对观测平面进行特殊剖分。因此, 通过详细推导, 剖析解析奇点存在原因, 并寻找新途径构建新的无解析奇点公式显得十分必要。

	Figure Option View Download New Window
	图1 矩形棱柱体与外部场点相互关系示意

1 问题的提出

在众多重力梯度张量分量正演数值模拟和反演研究的文献^{[1, 2, 3, 4, 5]}中, 如图1所示, 矩形棱柱体不含积分限的解析公式为

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρarctan \frac{(η - y) (ζ - z)}{(ξ - x) r}, \\ U_{xy} = G_{0} ρln [(ζ - z) + r], \\ U_{xz} = G_{0} ρln [(η - y) + r] \\ U_{yy} = - G_{0} ρarctan \frac{(ξ - x) (ζ - z)}{(η - y) r}, \\ U_{yz} = G_{0} ρln [(ξ - x) + r], \\ U_{zz} = - G_{0} ρarctan \frac{(ξ - x) (η - y)}{(ζ - z) r}, (1) \end{matrix} \end{matrix}$

其中, r= $\begin{matrix} \sqrt{{(ξ - x)}^{2} + {(η - y)}^{2} + {(ζ - z)}^{2}} \end{matrix}$ , G₀=6.672× 10^-¹¹ m³/(kg· s²), 为万有引力常数。

式(1)的离散计算公式为

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \arctan \frac{(η_{j} - y) (ζ_{k} - z)}{(ξ_{i} - x) r_{ijk}}, \\ U_{xy} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \ln [(ζ_{k} - z) + r_{ijk}], \\ U_{xz} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \ln [(η_{j} - y) + r_{ijk}], \\ U_{yy} = - G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \arctan \frac{(ξ_{i} - x) (ζ_{k} - z)}{(η_{j} - y) r_{ijk}}, \\ U_{yz} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \ln [(ξ_{i} - x) + r_{ijk}], \\ U_{zz} = - G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \arctan \frac{(ξ_{i} - x) (η_{j} - y)}{(ζ_{k} - z) r_{ijk}}, (2) \end{matrix} \end{matrix}$

其中:r_ijk= $\begin{matrix} \sqrt{(ξ_{i} {- x)}^{2} + (η_{j} {- y)}^{2} + (ζ_{k} {- z)}^{2}} \end{matrix}$ , s=(-1)^i+j+k。

公式(1)和公式(2)形式简洁、对称, 但在正演数值模拟中, 当观测平面网格剖分点与矩形棱柱体四个角点在观测平面投影重合时, U_xx和U_yy在该点无计算值, 为解析奇点, 如图2所示。在实际应用中, 为了避免奇点出现, 通常在观测平面采用特殊的剖分方法, 避开棱柱体角点在观测平面的投影点。

	Figure Option View Download New Window
	图2 U_xx和U_yy分量含奇点正演数值模拟结果

2 奇点存在原因分析

以下从基本公式出发推导矩形棱柱体解析公式, 以便说明奇点存在的原因。众所周知, 根据牛顿万有引力定律, 在笛卡尔坐标系中体密度为ρ (ξ , η , ζ )的任意三度体在其外部Q(x, y, z)点处的重力位(严格上为引力位)^[8]可表示为

$\begin{matrix} U (x, y, z) = G_{0} \underset{V}{\int \int \int} \frac{ρ (ξ, η, ζ)}{r} dξdηdζ, (3) \end{matrix}$

重力梯度是重力位的二阶导数, 写成张量形式为

$\begin{matrix} [\begin{matrix} U_{xx} & U_{xy} & U_{xz} \\ U_{yx} & U_{yy} & U_{yz} \\ U_{zx} & U_{zy} & U_{zz} \end{matrix}] = [\begin{matrix} \frac{\partial^{2} U}{\partial x^{2}} & \frac{\partial^{2} U}{\partial x \partial y} & \frac{\partial^{2} U}{\partial x \partial z} \\ \frac{\partial^{2} U}{\partial y \partial x} & \frac{\partial^{2} U}{\partial y^{2}} & \frac{\partial^{2} U}{\partial y \partial z} \\ \frac{\partial^{2} U}{\partial z \partial x} & \frac{\partial^{2} U}{\partial z \partial y} & \frac{\partial^{2} U}{\partial z^{2}} \end{matrix}] 。 (4) \end{matrix}$

对于图1所示密度均匀矩形棱柱体而言, 将式(3)代入式(4), 并令ρ (ξ , η , ζ )=ρ 、ξ -x=x'、η -y=y'、ζ -z=z', 则dξ =dx'、dη =dy'、dζ =dz', 可得常用的6个重力梯度张量分量计算公式:

$\begin{matrix} \begin{matrix} U_{xx} = G_{0} \underset{V}{\int \int \int} ρ \frac{3 x'^{2} - r^{2}}{r^{5}} dx'dy'dz', \\ U_{xy} = G_{0} \underset{V}{\int \int \int} ρ \frac{3 x'y'}{r^{5}} dx'dy'dz', \\ U_{xz} = G_{0} \underset{V}{\int \int \int} ρ \frac{3 x'z'}{r^{5}} dx'dy'dz', \\ U_{yy} = G_{0} \underset{V}{\int \int \int} ρ \frac{3 y'^{2} - r^{2}}{r^{5}} dx'dy'dz', \\ U_{yz} = G_{0} \underset{V}{\int \int \int} ρ \frac{3 y'z'}{r^{5}} dx'dy'dz', \\ U_{zz} = G_{0} \underset{V}{\int \int \int} ρ \frac{3 z'^{2} - r^{2}}{r^{5}} dx'dy'dz' 。 (5) \end{matrix} \end{matrix}$

根据重磁位场泊松公式, 式(5)与三分量磁异常的核函数完全一致^[9]。实际应用中, 仅需计算矩形棱柱体上半无源空间的重力梯度张量。同时, 由于U_xx、U_yy、U_zz及U_xy、U_xz、U_yz分别在形式上具有相似特征, 只有推导出其中任意两个, 然后再进行简单变量替换便可获得其他几个张量分量的积分结果。

对U_xx采用分步积分和变量替换方法先计算其三重不定积分:

$\begin{matrix} \begin{matrix} U_{xx} = G_{0} ρ \underset{V}{\int \int \int} \frac{3 x'^{2} - r^{2}}{r^{5}} dx'dy'dz' = \\ G_{0} ρ \underset{V}{\int \int \int} \frac{2 x'^{2} - y'^{2} - z'^{2}}{(x'^{2} + y'^{2} + z'^{2})^{\frac{5}{2}}} dx'dy'dz' 。 (6) \end{matrix} \end{matrix}$

详细推导过程如下:首先对x'进行积分, 在上半无源空间z'=ζ -z> 0, 利用变量替换, 令y'²+z'²=a²(a > 0), x'=a· tan t, (- $\begin{matrix} \frac{π}{2} \end{matrix}$ < t< $\begin{matrix} \frac{π}{2} \end{matrix}$ ), 则dx'=a· sec²t dt, 代入得

$\begin{matrix} \begin{matrix} U_{xx} = G_{0} ρ \iint [\int \frac{2 x'^{2} - y'^{2} - z'^{2}}{(x'^{2} + y'^{2} + z'^{2})^{\frac{5}{2}}} dx'] dy'dz' = \\ G_{0} ρ \iint [\int \frac{2 a^{2} ta n^{2} t - a^{2}}{(a^{2} ta n^{2} t + a^{2})^{\frac{5}{2}}} \cdot ase c^{2} tdt] dy'dz' = \\ G_{0} ρ \iint [\int \frac{1}{a^{2}} (3 si n^{2} t - 1) dsin t] dy'dz' = \\ G_{0} ρ \iint \frac{1}{a^{2}} (si n^{3} t - sin t) dy'dz', (7) \end{matrix} \end{matrix}$

又:x'=a· tan t, 得到sint= $\begin{matrix} \frac{x'}{\sqrt{x'^{2} + a^{2}}} \end{matrix}$ , 带入式(7)得:

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρ \iint \frac{x'}{(x'^{2} + y'_{2} + z'^{2})^{\frac{3}{2}}} dy'dz' = \\ - G_{0} ρx' \iint \frac{1}{(x'^{2} + y'_{2} + z'^{2})^{\frac{3}{2}}} dy'dz' 。 (8) \end{matrix} \end{matrix}$

同理, 令x'²+z'²=b²(b> 0), y'=b· tan t, (- $\begin{matrix} \frac{π}{2} \end{matrix}$ < t< $\begin{matrix} \frac{π}{2} \end{matrix}$ ), 则dy'=b· sec²t dt, 带入式(8)得

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρx' \int [\int \frac{1}{(x'^{2} + y'^{2} + z'^{2})^{\frac{3}{2}}} dy'] dz' = \\ - G_{0} ρx' \int [\int \frac{1}{b^{3} se c^{3} t} \cdot bse c^{2} tdt] dz' = \\ - G_{0} ρx' \int \frac{sin t}{b^{2}} dz', (9) \end{matrix} \end{matrix}$

又:y'=btan t, 得到sin t= $\begin{matrix} \frac{y'}{\sqrt{y'^{2} + b^{2}}} \end{matrix}$ , 代入可得

$\begin{matrix} U_{xx} = - G_{0} ρx'y' \int \frac{1}{(x'^{2} + z'^{2}) (x'^{2} + y'^{2} + z'^{2})^{\frac{1}{2}}} dz' 。 (10) \end{matrix}$

参照前两步积分方法, 利用变量替换, 令x'²+y'²=c², 当x', y'不同时为零时, c> 0, z'=ctan t, (- $\begin{matrix} \frac{π}{2} \end{matrix}$ < t< $\begin{matrix} \frac{π}{2} \end{matrix}$ ), 则dz'=c· sec²tdt, 代入式(10)得

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρx'y' \int \frac{1}{(x'^{2} + z'^{2}) (x'^{2} + y'^{2} + z'^{2})^{\frac{1}{2}}} dz' = \\ - G_{0} ρx'y' \int \frac{1}{(x'^{2} + c^{2} ta n^{2} t) csec t} \cdot cse c^{2} tdt = \\ - G_{0} ρx'y' \int \frac{cos t}{x'^{2} co s^{2} t + c^{2} si n^{2} t} dt = \\ - G_{0} ρx'y' \int \frac{1}{x'^{2} + y'^{2} si n^{2} t} dsin t = \\ - G_{0} ρ \int \frac{1}{1 + {(\frac{y'}{x'} sin t)}^{2}} d (\frac{y'}{x'} sin t) = \\ - G_{0} ρarctan \frac{y'sin t}{x'}, (11) \end{matrix} \end{matrix}$

又:z'=c· tan t, 得到sin t= $\begin{matrix} \frac{z'}{\sqrt{z'^{2} + c^{2}}} \end{matrix}$ , 带入式(11)可得

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρarctan \frac{y'z'}{x' \sqrt{x'^{2} + y'^{2} + z'^{2}}} = \\ - G_{0} ρarctan \frac{y'z'}{x'r} 。 (12) \end{matrix} \end{matrix}$

由于对角分量的公式相似, 简单交换积分次序并利用类似变量替换便可得到其他两个分量公式, 因此按照上述积分方法得到对角分量解析表达式为

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρarctan \frac{y'z'}{x'r}, \\ U_{yy} = - G_{0} ρarctan \frac{x'z'}{y'r}, \\ U_{zz} = - G_{0} ρarctan \frac{x'y'}{z'r} 。 (13) \end{matrix} \end{matrix}$

上述积分方法, 在推导过程中前两步积分都是严格成立的, 但第三步对z'积分过程并不严谨。首先, 在上半无源空间, 矩形棱柱体顶面角点在计算平面上的投影点存在x'、y'同时为零的情况, 即:x'²+y'²=c²=0, 此时, 不满足变量替换条件, 由此得出的积分结果无意义, 为“ 解析奇点” ; 其次, 在积分倒数第二步分子分母同除x'² 亦不严谨, 因为在矩形棱柱体顶面角点在计算平面上的投影点存在x'=0的情况。

3 无解析奇点公式

如上所述, 利用前述变量代换方法对z'的积分会产生解析奇点, 为此借鉴了郭志宏等人提出的长方体无解析奇点的积分方法^[7]:

$\begin{matrix} \begin{matrix} U_{xx} = - G_{0} ρ \int \frac{x'y'}{(x'^{2} + z'^{2}) (x'^{2} + y'^{2} + z'^{2})^{\frac{1}{2}}} dz' = \\ - G_{0} ρ \int \frac{x'y' \cdot (r + z'^{2}) / r}{(x'^{2} + z'^{2}) \cdot {(r + z')}^{2}} dz' = \\ G_{0} ρ \int \frac{d (\frac{x'y'}{x'^{2} + rz' + z'^{2}})}{1 + {(\frac{x'y'}{x'^{2} + rz' + z'^{2}})}^{2}} = \\ G_{0} ρarctan \frac{x'y'}{x'^{2} + rz' + z'^{2}}, (14) \end{matrix} \end{matrix}$

U_yy与U_zz积分公式具有相同形式, 同理可得

$\begin{matrix} U_{yy} = G_{0} ρarctan \frac{x'y'}{y'^{2} + rz' + z'^{2}} 。 (15) \end{matrix}$

虽然U_zz积分公式形式相同, 但采用该方法推导时, 前两步应该先对z'进行积分, 不存在最后对z'积分的情况。根据拉普拉斯方程:U_xx+U_yy+U_zz=0, 可得对角分量解析表达式:

$\begin{matrix} \begin{matrix} U_{zz} = - G_{0} ρ [\arctan \frac{x'y'}{x'^{2} + rz' + z'^{2}} + \\ \arctan \frac{x'y'}{y'^{2} + z'r + z'^{2}}], (16) \end{matrix} \end{matrix}$

由于: $\begin{matrix} \frac{x'y'}{x'^{2} + rz' + z'^{2}} \end{matrix}$ · $\begin{matrix} \frac{x'y'}{y'^{2} + rz' + z'^{2}} \end{matrix}$ < 1, 利用反三角函数性质:artan x+arctan y=arctan $\begin{matrix} \frac{x + y}{1 - xy} \end{matrix}$ , (xy< 1), 可得

$\begin{matrix} U_{zz} = - G_{0} ρarctan \frac{x'y'}{z'r} 。 (17) \end{matrix}$

综上可得, 重力梯度张量对角分量解析表达式为

$\begin{matrix} \begin{matrix} U_{xx} = G_{0} ρarctan \frac{x'y'}{x'^{2} + z'r + z'^{2}}, \\ U_{yy} = G_{0} ρarctan \frac{x'y'}{y'^{2} + z'r + z'^{2}}, \\ U_{zz} = - G_{0} ρar ctan \frac{x'y'}{z'r}, (18) \end{matrix} \end{matrix}$

式(18)在Bhattacharyya B K关于矩形棱柱体磁异常公式推导的文献^[10^-^11]中出现过类似形式。

重力梯度非对角张量分量U_xy、U_xz、U_yz分别在形式上也具有相似特征, 其积分过程以U_xy为例推导:

$\begin{matrix} \begin{matrix} U_{xy} = G_{0} ρ \underset{V}{\int \int \int} \frac{3 x'y'}{(x'^{2} + y'^{2} + z'^{2})^{\frac{5}{2}}} dx'dy'dz' = G_{0} ρ \cdot \\ \int \int [y' \int \frac{3}{{(\sqrt{x'^{2} + y'^{2} + z'^{2}})}^{4}} d \sqrt{x'^{2} + y'^{2} + z'^{2}}] dy'dz' = \\ - G_{0} ρ \iint \frac{y'}{{(\sqrt{x'^{2} + y'^{2} + z'^{2}})}^{3}} dy'dz' = \\ - G_{0} ρ \int [\int \frac{1}{{(\sqrt{x'^{2} + y'^{2} + z'^{2}})}^{2}} d \sqrt{x'^{2} + y'^{2} + z'^{2}}] dz' = \\ G_{0} ρ \int \frac{1}{\sqrt{x'^{2} + y'^{2} + z'^{2}}} dz' = \\ G_{0} ρln (z' + \sqrt{x'^{2} + y'^{2} + z'^{2}}) = \\ G_{0} ρln (z' + r), (19) \end{matrix} \end{matrix}$

对于U_xz、U_yz, 只需要简单地交换积分次序便可得到解析表达式, 因此, 矩形棱柱体重力梯度非对角张量分量解析表达式:

$\begin{matrix} \begin{matrix} U_{xy} = G_{0} ρln (z' + r), \\ U_{xz} = G_{0} ρln (y' + r), \\ U_{yz} = G_{0} ρln (x' + r) 。 (20) \end{matrix} \end{matrix}$

综合公式(19)和公式(21), 可得矩形棱柱体重力梯度张量分量无奇点解析公式:

$\begin{matrix} \begin{matrix} U_{xx} = G_{0} ρarctan \frac{(η - y) (ζ - z)}{{(ξ - x)}^{2} + (ζ - z) r + {(ζ - z)}^{2}}, \\ U_{xy} = G_{0} ρln [(ζ - z) + r], \\ U_{xz} = G_{0} ρln [(η - y) + r], \\ U_{yy} = G_{0} ρarctan \frac{(η - y) (ζ - z)}{{(η - y)}^{2} + (ζ - z) r + {(ζ - z)}^{2}}, \\ U_{yz} = G_{0} ρln [(ξ - x) + r], \\ U_{zz} = - G_{0} ρarctan \frac{(ξ - x) (η - y)}{(ζ - z) r}, (21) \end{matrix} \end{matrix}$

无解析奇点离散计算公式:

$\begin{matrix} \begin{matrix} U_{xx} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \\ \arctan \frac{(ξ_{i} - x) (η_{j} - y)}{(ξ_{i} {- x)}^{2} + (ζ_{k} - z) r_{ijk} + (ζ_{k} {- z)}^{2}}, \\ U_{xy} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \ln [(ζ_{k} - z) + r_{ijk}], \\ U_{xz} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \ln [(η_{j} - y) + r_{ijk}], \\ U_{yy} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \\ \arctan \frac{(ξ_{i} - x) (η_{j} - y)}{(η_{j} {- y)}^{2} + (ζ_{k} - z) r_{ijk} + (ζ_{k} {- z)}^{2}}, \\ U_{yz} = G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \ln [(ξ_{i} - x) + r_{ijk}], \\ U_{zz} = - G_{0} ρ \overset{2}{\sum_{i = 1}} \overset{2}{\sum_{j = 1}} \overset{2}{\sum_{k = 1}} s \cdot \arctan \frac{(ξ_{i} - x) (η_{j} - y)}{(ζ_{k} - z) r_{ijk}} 。 (22) \end{matrix} \end{matrix}$

4 模型试验及理论分析

为了检验理论公式的正确性, 选用Marshall M Rogers论文中的正演模型^[1], 其空间分布如图3所示, 坐标系z轴向下, 模型参数如下:体积50 m× 10 m× 6 m, x范围为100~150 m, y范围为120~130 m, z范围为47~53 m, 模型的中心坐标(125 m, 125 m, 50 m); 物性参数:密度差为1.5 g/cm³。

	Figure Option View Download New Window
	图3 正演模型空间分布示意

地面测网大小250 m× 250 m, 网格间距5 m× 5 m, x、y范围均为0~250 m, 节点数共51× 51。

	Figure Option View Download New Window
	图4 重力梯度张量异常正演数值模拟结果(无解析奇点正演结果)

图4、图5所示分别为新推导无奇点公式(22)和原有含奇点公式(2)计算结果, 对比可见, 后者除U_xx、U_yy分量含有4个奇点外, 二者在异常特征和异常强度上完全一致。Matlab正演数值模拟计算结果显示:精确到小数点后四位, 二者异常强度相同, 证明了新推导的无奇点公式的正确性。

	Figure Option View Download New Window
	图5 重力梯度张量异常正演数值模拟结果(含解析奇点正演结果)

5 结论

1)采用全新的推导方法详细推导了矩形棱柱体重力梯度张量异常含奇点正演计算公式, 结合推导过程分析了现有公式中解析奇点产生的原因。

2)借鉴文献[7]中关于长方体磁总场异常的推导方法, 推导了矩形棱柱体重力梯度张量异常无解析奇点计算公式, 通过坐标旋转, 该套公式适用于所有直立矩形棱柱体上方观测平面重力梯度张量异常的计算。

3)模型试验结果表明:重力梯度张量无奇点公式有效的消除了常用公式中的解析奇点, 与原有含奇点公式计算结果对比显示, 除奇点外二者计算结果一致(精确到小数点后四位), 证明了无奇点正演计算公式的正确性。

The authors have declared that no competing interests exist.

参考文献

View Option

[1]	Marshall M Rogers. An Investigation into the Feasibility of using a Modern Gravity Gradient Instrument for Passive Aircraft Navigation and Terrain Avoidance[C]. Thesis for Degree of Master, 2009. [本文引用:3]
[2]	陈召曦, 孟小红, 郭良辉, 等. 基于GPU并行的重力、重力梯度三维正演快速计算及反演策略[J]. 地球物理学报, 2012, 55(12): 4069-4077. [本文引用:2]
[3]	Mark Pilkington. Evaluating the utility of gravity gradient tensor components[J]. Geophysics, 2014, 79(1): 1-14. [本文引用:2]
[4]	陈闫, 李桐林, 范翠松, 等. 重力梯度全张量数据三维共轭梯度聚焦反演[J]. 地球物理学进展, 2014, 29(3): 1133-1142. [本文引用:2]
[5]	杨娇娇, 刘展, 陈晓红, 等. 基于深度加权的重力梯度张量数据的3D聚焦反演[J]. 地质世界, 2015, 34(1): 210-218. [本文引用:2]
[6]	Donald, Plouff. Gravity and magnetic fields of polygonal prisms and application to magnetic terrain corrections[J]. Geophysics, 1976, 41(4): 727-741. [本文引用:1]
[7]	郭志宏, 管志宁, 熊盛青. 长方体Δ T场及其梯度场无解析奇点理论表达式[J]. 地球物理学报, 2004, 47(6): 1131-1138. [本文引用:2]
[8]	曾华霖. 重力场与重力勘探[M]. 北京: 地质出版社, 2005. [本文引用:1]
[9]	Manik Talwani. Computation with the help of a digital computer of magnetic anomalies caused by bodies of arbitrary shape[J]. Geophysics, 1965, 30(5): 797-817. [本文引用:1]
[10]	Bhattacharyya B K. Magnetic anomalies due to prism-shaped bodies with arbitrary polarization[J]. Geophysics, 1964, 29(4): 517-531. [本文引用:1]
[11]	Bhattacharyya B K. A method for computing the total magnetization vector and the dimensions of a rectangular block-shaped body from magnetic anomalies[J]. Geophysics, 1966, 31(1): 74-96. [本文引用:1]
[12]	Dezso Nagy. The gravitational attraction of a right rectangular prism[J]. Geophysics, 1966, 31(2): 362-371. [本文引用:1]
[13]	骆遥, 姚长利. 长方体磁场及其梯度无解析奇点表达式理论研究[J]. 石油地球物理勘探, 2007, 42(6): 714-719. [本文引用:1]
[14]	骆遥, 姚长利, 薛典军, 等. 2. 5D地质体重磁异常无解析奇点正演计算研究[J]. 石油地球物理勘探, 2009, 44(4): 487-493. [本文引用:1]

2009

0.0

... 笔者研究发现,目前被广泛引用的矩形棱柱体重力梯度张量异常计算公式^[1,2,3,4,5]均源自Donald Plouff^[6],该公式虽然简洁、对称,但文献中仅给出了结果,没有详细地推导过程,且存在解析奇点^[7],在实际应用中需要对观测平面进行特殊剖分 ...