第４５卷　第３期２０２３年５月
物探化探计算技术
ＣＯＭＰＵＴＩＮＧＴＥＣＨＮＩＱＵＥＳＦＯＲＧＥＯＰＨＹＳＩＣＡＬＡＮＤＧＥＯＣＨＥＭＩＣＡＬＥＸＰＬＯＲＡＴＩＯＮ
Ｖｏｌ．４５　Ｎｏ．３
Ｍａｙ２
０２３　
收稿日期：２０２２ ０１ ２８第一作者：崔亚彤（１９９３－），女，博士，主要从事地球物理数据处理及相关研究工作，Ｅ ｍａｉｌ：ＹａｔｏｎｇＣｕｉ＠ｅｍａｉｌ．ｃｕｇ
ｂ．ｅｄｕ．ｃｎ。文章编号：１００１ １７４９（２０２３）０３ ０３５１ ０９球坐标系变密度界面正演及反演方法研究
崔亚彤１，王胜侯２
（１．天津市勘察设计院集团有限公司，天津　３００１９１；２．中国地质大学（武汉）　资源学院，
武汉　４３００７４）摘　要：密度界面反演方法可推断构造界面起伏形态，对于分析区域地质构造、地球深部构造和划分构造单元具有重要意义。当涉及到区域性乃至全球尺度问题时，
通常采用球坐标系的密度界面反演方法。传统球坐标系密度界面反演方法在地表观测面的计算精度较低，且未能考虑地下密度的纵横向变化。这里针对上述问题，开展球坐标系变密度界面正演及反演方法研究，通过理论推导和高斯－勒让德积分算法改进，给出球坐标系变密度界面高精度正演方法，运用迭代优化算法而给出球坐标系变密度界面反演方法。最后利用球壳模型和起伏界面数据试验验证该方法的有效性，适用于大区域、地表观测面的变密度界面反演。关键词：重力异常；球坐标系；变密度界面正演；变密度界面反演
中图分类号：Ｐ６３１．１  文献标志码：Ａ  犇犗犐：１０．３９６９／ｊ
．ｉｓｓｎ．１００１ １７４９．２０２３．０３．０９０　引言
密度界面反演是重磁定量解释的一个重要手段，
用于推断基底面、莫霍面、岩石圈－软流圈边界等构造界面形态，长期在国内、外应用较广。频率域密
度界面反演方法是一种常用方法，可以快速给出密度界面深度和起伏形态。Ｐａｒｋｅｒ－Ｏｌｄｅｎｂｕｒｇ方
法是目前最为经典、流行的频率域界面反演方法［１］
。
频率域的密度界面反演方法虽然计算简单、快速，但该方法基于笛卡尔坐标系，且仅适用于小区域尺度的反演。当涉及到区域性乃至全球尺度问题，采用球坐标系密度界面反演方法，可以缓解大区域尺度
的地球曲率问题［２－６
］，对于区域性乃至全球尺度的
深部构造界面反演、区域构造地质解释等，具有重要的意义。针对球坐标系密度界面反演问题，通常是
基于直接迭代法［２－３］或正则化非线性重力反演方
法［４
］，前者无需进行非线性矩阵方程的运算，计算相
对简单、高效。
直接迭代法根据研究区实测重力异常和界面上下密度差，通过不断对估计密度界面深度进行正演计算，以达到理论异常和实测异常的最小化，最终可得到研究区密度分界面深度起伏情况。因此高精度的球坐标系正演算法，对于密度界面反演起着至关重要的作用。通过利用球坐标系与笛卡尔坐标系的转换关系，可以将常规基于笛卡尔坐标系的直立长方体正演公式转换到球坐标系中，原来的直立长方体就变为球坐标系中的球面棱柱体，称之为Ｔｅｓｓｅ ｒｏｉｄ单元体。对于整个研究区来说，
观测面上任意一点的重力异常值，即是由正演计算得到的地下Ｔｅｓｓｅｒｏｉｄ单元体在该点处引起的重力异常的累加求和，遍历整个观测网即可得到由密度界面引起的重力异常，这即是球坐标系密度界面正演方法。
球坐标系中的均质Ｔｅｓｓｅｒｏｉｄ单元体正演积分
公式不具有解析解，为此很多学者利用数值积分方
法得到其数值解，常用的数值积分方法主要是三维高斯－勒让德数值积分（３ＤＬｅｇ
ｅｎｄｒｅ－Ｇａｕｓｓｑｕａｄｒａｔｕｒｅ，３Ｄ－ＧＬＱ）方法［７－８
］和泰勒公式（Ｔａｙ
ｌｏｒ）数值积分方法［９］
。通过引入自适应网格剖分方法［
１０］
，虽然可以在一定程度上提高正演精度，但大多数算法仍然只适用于高空（卫星）数据，不适用于地表观测面。
在密度界面正演和反演中，密度参数是关乎算法准确性的最主要参数之一。在实际应用中，大区域尺度研究区的地层密度在水平和铅垂方向的变化均很复杂，利用常数作为密度差异进行反演的结果显然不够准确。因此，
在界面反演过程中引入密度变化模型，可以更好地模拟实际地层密度分布情
况［１１－１２］
，进一步提高反演精度。球坐标系的变密度界面反演方法成果较少，Ｌｉｎ等［８］推导出了线性变
密度Ｔｅｓｓｅｒｏｉｄ单元体的正演公式，但该正演公式是利用Ｔａｙｌｏｒ数值积分方法得到的数
值解，Ｔａｙｌｏｒ的高阶展开公式虽然计算精度高，但公式繁琐，计算量大，利用低阶展开公式又难以满足正演精度。
这里针对大区域尺度变密度界面正演及反演问题展开理论推导和技术改进，
研究给出球坐标系变密度界面正演和反演方法，为区域构造界面研究提供技术支撑。
１　球坐标系变密度界面反演方法
球坐标系密度界面反演是推断大区域基底面、莫霍面、岩石圈－软流圈边界等构造界面形态的重要方法。界面反演主要通过反复的正演计算及迭代拟合来获取界面起伏情况。因此高精度的正演算法，对反演计算起着至关重要的作用。
１．１　球坐标系变密度界面正演方法
球坐标系（图１（ａ））将三个方向矢量分别定义为经度方向、
纬度方向、以及沿半径指向球心的径向。基于球坐标系的重力异常正演方法，即是将测区（计算区域）以经度间距Δλ，纬度间距Δφ的两组直线划分成一系列网格，每个结点就代表一个横截面为Δλ×Δφ、
深度为Δ狉犻犼的Ｔｅｓｓｅｒｏｉｄ单元体（图１）
，图１　球坐标系剖分网格示意图
Ｆｉｇ．１　Ｇｒｉｄｄｉｖｉｓｉｏｎｕｎｄｅｒｓｐｈｅｒｉｃａｌｃｏｏｒｄｉｎａｔｅｓｙ
ｓｔｅｍ各单元体的密度可以相同，也可不同，那么每个测点处的重力异常即是所有Ｔｅｓｓｅｒｏｉｄ单元体在该点处重力异常的叠加。
Ｈｅｃｋ等［９］
根据笛卡尔坐标系与球坐标系的坐
标关系，
将笛卡尔坐标系中的重力异常公式转换到球坐标系中，则基于球坐标系的Ｔｅｓｓｅｒｏｉｄ单元体
重力异常积分公式可表示为式（１）。Δ犵（狉，φ，λ）＝Ｇ∫λ２λ
１
∫φ２φ１
∫
狉２
狉１
ρ狉＇２
（狉－狉＇ｃｏｓψ）ｃｏｓφ＇犾３ｄ狉＇ｄφ＇ｄλ＇（１
）其中：犾为观测点（狉，φ，λ）与源内点（狉＇，φ＇，λ＇）之间的距离；ψ为观测点和源内点的位置向量之间的夹
角［
９
］犾＝狉２＋狉＇２
－２狉狉＇ｃｏｓ槡ψ
ｃｏｓψ＝ｓｉｎφｓｉｎφ＇＋ｃｏｓφｃｏｓφ
＇ｃｏｓ（λ＇－λ）（２
）其中：Ｇ为万有引力常数，Ｇ＝６．６７×１０－１１ｍ３
／ｋｇ·ｓ２；ρ为密度差（ｋｇ／ｍ３
）；（狉１，狉２）、（φ１，φ２）、（λ１，λ２）
分别为Ｔｅｓｓｅｒｏｉｄ单元体径向、纬度方向、经度方向积分上下限，且有Δ狉＝狉２－狉１，Δφ＝φ２－φ１，Δλ＝λ２－λ１。
在处理实际问题中，地下密度通常是随空间而变化的。为了模拟地下实际密度分布情况，前人基于笛卡尔坐标系提出多种密度－深度函数模型［１１－１２］。在笛卡尔坐标系中，通常称深度方向密度变化为纵向变密度，为了加以区分，将其称之为径向变密度。假设Ｔｅｓｓｅｒｏｉｄ单元体密度随径向线性变化时，即公式（１）中ρ（狉＇）＝ρ０＋犪狉＇。因此径向线性变密度重力异常正演公式可表示为式（３
）。Δ犵（狉，φ，
λ）＝Ｇ∫λ２λ１∫φ２φ１
∫
狉２
狉１ρ（狉＇）狉＇２（狉－狉＇ｃｏｓψ）ｃｏｓφ＇犾３ｄ狉＇ｄφ＇ｄλ＇＝Ｇ∫λ２λ１
∫
φ２
φ１
ｃｏｓφ＇ｄφ
＇ｄλ＇·∫
狉２
狉１
ρ０狉＇２
（狉－狉＇ｃｏｓψ）ｄ狉＇犾３
＋∫
狉２
狉１
犪狉＇３
（狉－狉＇ｃｏｓψ）ｄ狉＇犾熿
燀燄燅
３
（３
）２
５３    物探化探计算技术４５卷
其中：犪、ρ
０为线性参数，根据三维密度模型可直接获取。根据Ｇｒａｄｓｈｔｅｙｎ和Ｒｙ
ｚｈｉｋ整理的积分公式［１３
］，将式（３）
对狉进行积分，得到狉方向的解析解，进一步推导为式（４
）。Δ犵（狉，φ，λ）＝犌ρ０
狉
∫λ２λ１∫
φ２
φ１ｃｏｓφ＇ｄφ＇ｄλ＇狉＇３
犾－犾（狉＇＋３狉ｃｏｓψ）－狉２（３ｃｏｓ２
ψ－１）×ｌｎ狘２狉＇－２狉ｃｏｓψ＋２犾｛｝
狘狉２狉１
＋
Ｇ犪
∫λ２λ１
∫
φ２φ１
ｃｏｓφ＇ｄφ＇ｄλ＇１
犾（－１２ｃｏｓψ）狉＇３＋（狉－５２狉ｃｏｓ２ψ）狉＇２＋（１５狉２ｃｏｓ３ψ－１３２狉２ｃｏｓψ）狉＇＋（２狉３－１５２狉３ｃｏｓ２
ψ熿燀燄燅）＋３２狉２ｃｏｓψ（３－５ｃｏｓ２
ψ）ｌｎ狘２狉＇－２狉ｃｏｓψ＋２犾烅烄烆烍烌烎狘狉２狉
１
（４
）  再利用二维高斯－勒让德（２ＤＬｅｇｅｎｄｒｅ－Ｇａｕｓｓｑｕａｄｒａｔｕｒｅ
，
２Ｄ－ＧＬＱ）积分近似，式（４）解析解可近似为数值解，基于球坐标系的径向线性变密度界面正演公式为式（５
）。Δ犵（狉，φ，λ）＝Ｇρ０狉ΔλΔφ４∑犖
φ犼＝１∑犖λ
犻＝１
狑φ犼狑λ犻ｃｏｓ槇φ犼
狉＇３槇犾犻犼－槇犾犻犼（狉＇＋３狉ｃｏｓ槇ψ犻犼）－狉２（３ｃｏｓ２槇ψ犻犼－１）×ｌｎ狘２狉＇－２狉ｃｏｓ槇ψ犻犼＋２槇犾犻犼烅烄烆烍烌烎狘狉２狉１
＋Ｇ犪ΔλΔφ４∑犖
φ
正则化坐标
犼＝１∑犖λ
犻＝１狑φ犼狑λ犻ｃｏｓ槇φ犼１槇犾犻犼
（－１２ｃｏｓ槇ψ犻犼）狉＇３＋（狉－５２
狉ｃｏｓ２槇ψ犻犼）狉＇２＋（１５狉２ｃｏｓ３槇ψ犻犼－１３２狉２ｃｏｓ槇ψ犻犼）狉＇＋（２狉３
－１５２
狉３ｃｏｓ２槇ψ犻犼熿燀燄燅
）＋
３２狉２ｃｏｓ槇ψ犻犼（３－５ｃｏｓ２槇ψ犻犼）ｌｎ狘２狉＇－２狉ｃｏｓ槇ψ犻犼＋２槇犾犻犼烅烄烆烍
烌烎狘狉２狉
１
（５
）其中，槇
犾犻犼＝狉２＋狉＇２－２狉狉＇２ｃｏｓ槇ψ
犻槡犼ｃｏｓ槇ψ犻犼＝ｓｉｎφｓｉｎ槇φ犼＋ｃｏｓφｃｏｓ槇φ
犼ｃｏｓ（槇
λ犻－λ）槇
λ犻＝（λ２－λ１２）λ犻＋（λ２＋λ１２
）槇φ犼＝（φ２－φ１２）φ犼＋（φ２＋φ１２烅烄烆
）狑λ犻＝２（１－λ２犻）［犘狀＇（λ犻）
］２
狑φ
犼＝２（１－φ２犼）［犘狀＇（φ
犼）］烅
烄烆２
（犖ψ，犖λ）为２Ｄ－ＧＬＱ积分节点个数；（槇ψ
犼，槇
λ犻）为球坐标系中高斯节点的坐标信息；（狑φ犼，狑λ犻）为高斯－
勒让德系数；犘狀（
·）为勒让德多项式。为了进一步提高计算精度，使该方法适用于地表观测面，可以从两个思路解决这个问题，一是增加ＧＬＱ的高斯节点数，根据ＧＬＱ积分原理，
当高斯节点（犖ψ，犖λ）
越多时，数值解就越接近解析解，但是计算时间呈指数上升，
需要付出的计算代价巨大。另一个提高精度的方法就是自适应网格剖分方
法［８，１０
］，即当计算点接近源点时，将Ｔｅｓｓｅｒｏｉｄ单元
体再细剖分为几个小
的单元体（图２），然后在累加求和作为这个大单元体的重力异常。而当计算点远离源点一定距离时，则不进行剖分。
图２　自适应网格剖分示意图
（修改自Ｕｉｅｄａｅｔａｌ．［１
０］
）Ｆｉｇ．２　Ｓｃｈｅｍａｔｉｃｄｉａｇｒａｍｏｆａｄａｐ
ｔｉｖｅｇｒｉｄｄｉｖｉｓｉｏｎ（ｍｏｄｉｆｉｅｄｆｒｏｍＵｉｅｄａｅｔａｌ．［１０
］）
与前人提出的剖分方法［８］
不同的是，本文剖分
不考虑径向剖分，因为径向可求得解析解，所以仅考虑纬度和经度方向的剖分。通过下述条件式，可以判断网格是否进行剖分：
犾φ＝狉２（φ
２－φ１）犾λ＝狉２ｃｏｓφ
１（λ２－λ１）（６
）犾犾犱
≥犇　犱∈（φ，λ）（７）其中：犾φ、犾
λ为纬度方向和经度方向弧距；Ｄ为正数，称作距离尺寸比［１０］。笔者在利用自适应网格剖分
进行剖分时，发现存在无限次递归情况，致使程序卡顿。因此，为了控制剖分密集度，防止计算冗余，给出一个Ｆｌａｇ阈值来控制剖分次数，
一般设置为５。３５３３期
崔亚彤，等：球坐标系变密度界面正演及反演方法研究
自适应网格剖分方法步骤为：①判断犾与犾φ、犾
λ是否满足公式（６），如满足则无需进行网格剖分，如不满足则进行细致剖分；②若犾／犾λ＜Ｄ，则在经度方向进行剖分，若犾／犾φ＜Ｄ则在纬度方向进行剖分，若犾／犾λ＜Ｄ且犾／犾φ＜Ｄ
，则经纬度都需剖分；③经剖分迭代后，直至Ｆｌａｇ＝５或满足公式（７）。在正演过程中，除了需要引入自适应网格剖分方法，还需要在正演之前进行扩边，降低空间域计算过程中的边缘效应的影响。
据此，推导了球坐标系径向线性变密度Ｔｅｓｓｅ ｒｏｉｄ单元体正演公式，
在正演之前先进行扩边，为了适用于地表观测面，
进一步提高正演精度，在正演过程中需要采用自适应网格剖分方法进行网格剖分。对于变密度界面的重力异常正演，同样是将所有Ｔｅｓｓｅｒｏｉｄ单元体在任意观测点处重力异常进行累
加求和，最终可得到该密度界面所引起的重力异常。
若每个测点处，密度沿径向变化都一致时，那么使用线性拟合方法获得的每个测点处的犪、ρ０值都相等。但在实际数据处理中，先验密度模型在三个维度上均是变化的，在获得每个测点处的线性参数的同时，犪、ρ０值在横向上都是不相等的，因此在实现径向线性变密度的同时，
实现了横向变密度。然后再把犪、ρ
０值带入公式（５）计算即可实现三个方向的变密度界面正演。因此本文公式虽然是基于线性假设的Ｔｅｓｓｅｒｏｉｄ单元体的径向线性变密度正演公式，
但在处理实际问题时，
也可以同时实现横向变密度。１．２　球坐标系变密度界面反演方法
传统球坐标系重力反演方法通常适用于高空（卫星）观测面的重力数据，对地表观测面的并不适用，但在大区域尺度的研究区，地层密度通常存在横向和垂向的复杂变化，
宜采用变密度模型进行球坐标系密度界面的正演和反演。本文径向变密度界面
反演方法通过对Ｂｏｔｔ迭代方法进行优化［１
４］，通过不断地正演拟合，使得观测异常与计算理论异常达到误差允许范围内，
进而获得最终的界面起伏形态。针对径向变密度界面反演方法，搜集研究区内重力异常观测数据Δ犵０以及先验密度模型，利用该三维密度模型，通过使用线性拟合的方式获得每个测点处的犪、ρ０值。然后根据界面参考深度犺ｒｅｆ，结合线性公式ρ（狉）＝ρ
０＋犪（犚－犺ｒｅｆ），计算出每个测点处的密度值ρ，进而得到初始界面深度为式（８）。犺０＝
Δ犵０
２πＧρ
  狉０＝犚－犺０
（８）其中：狉０为初始界面所对应的径向深度；
犚
为地球半图３　球坐标系变密度界面迭代优化反演流程Ｆｉｇ．３　Ｆｌｏｗｏｆｖａｒｉａｂｌｅｄｅｎｓｉｔｙ
ｉｎｔｅｒｆａｃｅｉｔｅｒａｔｉｖｅ  ｏｐｔｉｍｉｚａｔｉｏｎｉｎｖｅｒｓｉｏｎｕｎｄｅｒ  ｓｐｈｅｒｉｃａｌｃｏｏｒｄｉｎａｔｅｓｙ
ｓｔｅｍ径，犚＝６３７１０００犿。将式（８）带入公式（５）进行正演计算，即可得到地表观测面的起伏界面的理论重力异常Δ犵１。利用Δ犵１与Δ犵０的残差更新初始界
面深度，即
犺１＝犺０＋
Δ犵０－Δ犵１
２πＧρ
，  狉１＝犚－犺１
（９）然后采用滑动窗口平滑的方法，将更新后的界
面起伏深度进行平滑处理。窗口平滑方法即对一定窗口内所有点的深度求取平均值，作为该窗口中心点的平滑深度值，窗口直径越大越平滑。利用平滑后界面深度再次进行正演计算，直至残差在误差允许范围之内，
输出最终的界面起伏深度。球坐标系变密度界面反演方法迭代流程如图３所示。
２　球壳模型数据正演试验
为了验证球坐标系线性变密度界面正演公式的正确性，采用球壳模型进行对比（图４）。假设其密度为随径向线性变化的，ρ（狉＇）＝ρ０＋犪狉＇，则密度线性变化的球壳重力异常表示为式（１０
）。Δ犵＝４πＧ狉２∫
狉２狉１
ρ（
狉＇）狉＇２
ｄ狉＇４πＧρ０３狉２（狉３２－狉３１）＋πＧ犪狉
２（狉４２－狉４１）（１０）假设球壳（图４）厚度为Δ狉＝４０ｋｍ，通过图４（ｂ）线性关系可以获得ρ０＝
６１５１，犪＝－０．００１，Ｔｅｓ ４
５３    物探化探计算技术４５卷
图４　径向线性变密度球壳模型
Ｆｉｇ．４　Ｒａｄｉａｌｌｉｎｅａｒｖａｒｉａｂｌｅｄｅｎｓｉｔｙｓｐ
ｈｅｒｉｃａｌｓｈｅｌｌｍｏｄｅｌ表１　径向线性变密度球壳模型重力异常对比Ｔａｂ．１　Ｇｒａｖｉｔｙａｎｏｍａｌｉｅｓｃｏｍｐ
ａｒｉｓｏｎｏｆｒａｄｉａｌｌｉｎｅａｒ  ｖａｒｉａｂｌｅｄｅｎｓｉｔｙｓｐ
ｈｅｒｉｃａｌｓｈｅｌｌ球壳公式２Ｄ－ＧＬＱ
２Ｄ－ＧＬＱ（
剖分）异常／ｍＧａｌ－６６６．４８－６５１．４０－６６６．４７百分误差／％
２．２６
０．０００６１
表２　不同高斯节点个数的计算精度效率对比Ｔａｂ．２　Ｃｏｍｐａｒｉｓｏｎｏｆｃａｌｃｕｌａｔｉｏｎａｃｃｕｒａｃｙａｎｄｅｆｆｉｃｉｅｎｃｙ
  ｗｉｔｈｄｉｆｆｅｒｅｎｔｎｕｍｂｅｒｏｆＧａｕｓｓｉａｎｎｏｄｅｓ高斯节点个数１２３４百分误差／％１．４０８０．２０６０．０５４０．００１计算时间／ｓ０．１３５０．２８３０．４７２０．８０８高斯节点个数５６７８百分误差／％０．０２５０．０３９０．０４７０．０５３计算时间／ｓ
１．３２８
２．１５３
３．４５４
５．６１１
ｓｅｒｏｉｄ网格间距１°×１°，网格大小３６１×１８１，
观测点位于地表，坐标为（６３７１０００，０，０），即犾＝０ｍ。分别采用变密度２Ｄ－ＧＬＱ、基于自适应网格剖分的变密度２Ｄ－ＧＬＱ以及球壳公式（式（１０））进行正演，高斯节点数犖φ＝犖λ＝４，正演结果如表１所示。结果显示，利用本文基于自适应网格剖分的变密度２Ｄ－ＧＬＱ积分方法计算得到的正演结果，
与球壳真实重力异常值基本相等，百分误差基本接近“０
”。同时，针对不同高斯节点个数进行正演精度以及计算时间对比试验（表２，图５），结果表明基于自适应剖分算法的正演计算效率，
随着高斯节点数的增加呈现指数上升的趋势，计算精度在高斯节点数为４时，精度最高，百分误差为０．００１％。
利用变密度球壳模型对球坐标系的变密度界面
图５　不同高斯节点个数的计算精度效率对比Ｆｉｇ．５　Ｃｏｍｐａｒｉｓｏｎｏｆｃａｌｃｕｌａｔｉｏｎａｃｃｕｒａｃｙａ
ｎｄ  ｅｆｆｉｃｉｅｎｃｙｗ
ｉｔｈｄｉｆｆｅｒｅｎｔｎｕｍｂｅｒ  ｏ
ｆＧａｕｓｓｉａｎｎｏｄｅｓ正演方法进行正演精度试验、评价。结果表明，本文球坐标系变密度界面正演方法计算精度高，适用于大区域、
地表观测面的密度界面的正演。３　密度界面模型数据反演试验
假设界面起伏模型（图６（ａ）），网格大小为１２５×７３，经度纬度方向网格间距为０．５°×０．５°，观测面位于地表，界面平均深度约３５ｋｍ。对于变密度正演计算，根据先验三维密度模型（图６（ｂ）），通过各个测点处的线性拟合，获取线性参数ρ０、犪，（图６（ｃ）和图６（ｄ）），通过图６（ｃ）、图６（ｄ）可以非常直观地看出，每个点处的线性参数ρ０、犪值均不同，因此可在实现径向线性变密度的同时，实现了横向变密度。
根据球坐标径向变密度正演公式（式（５）），上述界面起伏模型引起的重力异常可由公式（５）正演计算得到，图７（ａ）为变密度界面正演结果，界面深度深时，异常绝对值大，而界面深度浅时，异常
的绝对值相对小。
利用上述异常，分别对其进行常密度界面反演
５５３３期
崔亚彤，等：球坐标系变密度界面正演及反演方法研究

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。