剖分算法建模分析

1、剖分算法建模分析     引言地球重力场信息在地球构造研究、国土资源调查、矿产勘探以及国民经济等领域中发挥着重大作用重力的正问题与反问题是解释重力资料的核心研究内容，其中前者主要研究不同形状、产状和场源密度等场源体或地质体所引起的重力异常特征及分布等，其计算过程为正演计算；后者主要根据重力异常的分布来计算地质体或场源体的密度分布，其计算过程为反演计算目标地质体的重力建模建立了地质体的剩余密度值与地面重力异常值的关系，既可利用重力正演公式计算地质体所引起的重力异常，也可以应用于求解重力反问题因此，重力建模是重力正反演问题的关键环节，建模方法的优劣则直接影响

2、到重力正反演的计算精度随着对重磁资料解释的深入研究、重力测量精度的提高及计算机的快速发展，重力正反演研究由过去的二维逐步发展为三维对于三维连续密度变化的目标体，目前最常用的建模方法是块体方法，即将目标地质体剖分为规则的构造单元如长方体，计算所有单元产生的异常并叠加便得到目标地质体所产生的异常值楼海等给出了基于矩形网格模型的三维连续密度分布目标地质体的重力解析公式，并在实际应用中得到较好的结果，但当地质形状比较复杂时，其计算结果并不理想，如果加密构造单元，则又会增加计算量及存储空间盛国平采用水平截面法对非规则几何形状的目标体进行分割，并假设每一个质面层的密度值为常数，给出了较严密的数学解析表达式

3、张岭等针对截面为任意形状的二度体问题，利用二维剖分方法，将截面分割成若干三角形，二度体被分解成若干三棱柱的组合来计算二度体的重力异常值，这种方法表达了二维的连续密度变化，但仅限制于二度体的重力计算问题对于任意形状和变密度的三度体，现有重力建模方法均存在不足，为此，本文提出基于剖分算法的重力建模方法采用剖分算法将三维目标地质体分解为若干变密度四面体体元，推导基于四面体体元的重力正演公式；比较分析常规块体方法和本文方法应用于重力正演的计算效果，并采用共轭梯度法加密度约束条件对非规则形状变密度的地质体进行反演计算，验证本文重力建模方法的正确性和有效性剖分的重力建模地质体的建模地质体的建模步骤为：首先

4、根据地质体的三维坐标信息及物性信息等相关资料将目标地质体离散化为空间离散点；其次，采用合理的建模算法组织空间离散点，比如传统重力建模时采用的块体方法，地学领域、及剖分方法等；最后将所有体元结构（例如四面体、三棱柱体、六面体和长方体等）组合起来，根据体元顶点属性特征采用颜色或纹理对模型进行可视化对于规则形状的地质体建模，采用传统块体方法即可对地质体进行快速建模，并真实逼近模型结构及地质体的物性特征，但是对于非规则形状变密度地质体的建模，如褶皱、断层等，传统建模方法则会改变地质体的结构、密度分布和其他物性特征，所以本文引进剖分方法来逼近真实的地质体剖分算法由二维三角剖分算法演化而来，年为了限定二维

5、平面离散点的有效作用范围，首先定义了二维平面上的图，年由将图演化出了更易于分析应用的三角网图和三角网为目前普遍接受和广泛采用的分析研究区域离散数据的有力工具，已应用到石油勘探、地质、矿业、城市规划和环境监测等领域实现了网生成算法，对边界恢复算法进行改进，陈晓勇、李清泉等研究了四面体格网结构（）模型的生成算法剖分算法的剖分结果为四面体体元，四个顶点为空间离散点，包含坐标信息和物性特征四面体格网数据模型是二维三角形网（，）数据结构在三维空间上的扩展，与之相应的数据结构为四面体格网结构（）模型以三角剖分为基础，将目标空间用紧密排列但不重叠的非规则四面体的组合来表示四面体应具有如下特点：生成的四面体之

6、间无重叠部分；所形成的四面体组合可以覆盖整个目标体的三维空间；四面体的外接球不包含空间离散点集中除该四面体四个顶点外的任一点这些性质保证了组成四面体的三角形近似等边或等角，及四面体体元的组合更逼近真实的目标地质体为了比较常规块体方法及剖分算法对重力建模的影响，本文设计了两个目标地质体目标体为连续密度变化的矩形体，上顶面埋深为，范围为××（长×宽×高），八个角点的剩余密度值依次为、和，如图所示；目标体为倾斜台阶组合模型体，八个角点的剩余密度值与目标体相同，具体参数参见图对于目标体的建模，首先将目标体离散化为若干空间点位，再分别采用常规块体剖分方法与四面体剖

7、分方法对这些离散点位建模在、和方向取分块数为××，空间点数为××，对目标体的剖分结果见图比较图和图可知：两种方法的建模结果均能准确表达地质体，理论上得到的重力正演计算结果应该相差不大采用传统块体方法对目标体的建模结果见图当采用剖分算法对目标体建模时，如果直接利用空间离散点进行四面体构网，则会出现“跨越”现象，即左右倾斜台阶会出现错位结构，所以本文首先采用算法建立左右倾斜台阶接触面的数据结构，作为组合模型的特征面，然后进行整体四面体构网，其建模结果如图所示从图可以看出，具有约束的算法可以避免“跨越”现象，从而得到合理的三维模型比较图和图可知：块体方法的建模

8、结果与真实地质体有一定的差异，而剖分算法的建模结果能更好地逼近真实地质体，因此理论上两种方法的重力正演计算结果应该存在较大差异；块体方法的建模结果与分块数密切相关，而剖分算法受分块数的影响较小，随着图分块数的增加，块体方法的建模结果将逐渐逼近目标体的形状及密度分布，但计算量也将大幅度增加重力计算模型四面体体元为剖分算法的基本结构单元，非规则形状变密度地质体对地面点产生的重力异常为所有四面体体元产生的重力异常之和任意一个四面体体元（见图）对地面点产生的重力异常为式中，表示在四面体内任意点（，）的剩余密度，由四个顶点的剩余密度值线性插值得到，满足关系式：，（）式中，、和表示四面体体元四个顶点的剩余

9、密度值对于任意形状的四面体，利用式（）计算地面点的重力异常值，不易得到、和三个方向的积分域，因此也不利于重力异常的计算本文借鉴文献提出的水平截面法，首先将四面体切分为若干个水平面层，计算每一层所产生的重力异常，然后叠加得到一个四面体所引起的重力异常值与文献的不同之处在于：四面体四个顶点的密度是不同的，其每一个面层也是变密度的图为变密度的水平截面被划分为若干个三角形的示意图如图所示，为地面计算点在目标地质体的一个面层上的投影，则面层被分为个三角形，以为例，计算一个三角形对点产生的（）当线段平行于轴（）时，与的组成项变为，槡槡；当线段平行于轴（）时，与的组成项变为槡槡，（）、和三点在一条直线上，即

10、，则；（）根据格林公式求解与时，定义逆时针为正，要求、与、一致，即或时表示顺时针，系数，采用公式（）、（）、（）和（）计算出系数阵后，即可计算一个三角形质面对地面点产生的重力异常，按逆时针方向计算组成一个截面的个三角形对点产生的异常，将其叠加即为水平截面对地面点产生的重力异常假设四面体被分为了层水平截面，、和分别为相应埋深、和的水平截面所产生的重力异常值，则一个四面体对地面点产生的重力异常为综上所述，最终可建立目标体个剩余密度值与个地面重力异常的关系式，即算例分析目标体的正演计算分别采用常规块体剖分方法和剖分方法对目标体建模，并计算了位于目标体正上方范围为××的空间格网点的

11、重力异常，常规块体方法的正演计算公式可以参见文献，两种方法的计算结果如图所示从图可以看出：两种方法的计算结果相差较小，最大差值为×对目标体设计分块数为时，分别采用常规块体剖分方法和剖分方法对其建模，正演计算了位于目标体正上方处个格网点的重力异常，其差值的比较结果如表所示从表可以看出：随着分块数的增加，两种方法的计算结果并无明显差异，最大差值在×以内，均方误差在×以内因此，对于规则形状且连续密度变化的目标体，采用常规块体剖分方法和剖分方法均能获得较好的重力正演结果目标体的正演计算分别采用常规块体剖分方法和剖分方法对目标体建模，计算了与节相同区域的重力异常，计算结果如

12、图所示从图可以看出：对于非规则形状变密度的目标体，两种方法的计算结果差异较大，其差值在（）×范围之内，主要差异体现在双倾斜台阶组合模型体的公共斜平面上，其原因是常规块体方法改变了目标体的几何结构和物性特征，对地质体的建模不够准确类似于目标体，对目标体设计分块数为，分别采用常规块体方法和剖分方法对目标体建模并计算其产生的重力异常，两种方法的比较结果见表结果表明：随着分块数的增加，常规块体剖分方法对目标体的建模逐渐逼近真实地质体，使得两者计算结果的差值逐渐减小从表和表可以看出：平均值与均方误差量级相当，说明常规块体方法和剖分方法之间存在明显的系统误差，并且系统误差随着分块数的增加而减小这

13、主要是目标体的建模误差和式（）中划分层数的取值所导致根据文献和，的取值应根据研究区域及其地质条件来确定，本文的取值为目标体的反演计算剖分算法不仅可以应用于地质体的重力正演计算，也可以应用于地质体的三维物性反演针对地球物理的反演问题，和最早对连续线性问题做了深入研究；和引入深度加权函数，克服了三维反演结果的“趋附”现象；等将光滑度引入到正则化反演；在重磁反演求解大型线性方程组时，为了不降低计算效率及消耗过多计算机资源，采用共轭梯度法进行三维反演本文基于剖分算法，采用共轭梯度法加密度约束条件对目标体进行三维物性反演地质体剩余密度值与地面重力异常值之间的线性关系式（）往往是欠定的，不易得到真实物性反

14、演结果，为此本文增加了密度约束条件定义目标函数的表达式为这里为观测值，为实例模型结点剩余密度值，为该线性方程组的系数阵，、与分别对应于式（）的、与，和分别表示密度约束的下界和上界，在本算例的数值分别为采用最小二乘法求解式（）时，其矩阵求逆时间较长，而且有可能出现病态方程，因此本文利用共轭梯度法求解假定该区域的地质条件已知，背景密度为，以节基于本文方法正演计算的重力异常作为地面观测值，总共个地面重力异常值，其等值线图见图；采用个地面重力异常值反演目标体的个剩余密度值，其反演结果见图；利用密度反演结果计算得到的重力异常见图；图为图和图的重力异常之差，用于评价三维密度反演效果图的结果表明，采用共轭梯

15、度法加密度约束条件得到的反演结果非常接近于真实目标地质体的密度分布；从图可以看出：利用反演密度值计算的重力异常（见图）与模拟重力异常（见图相差很小，其量级为以上结果初步验证了本文的重力建模方法和三维密度反演方法的正确性和有效性结论与展望针对任意形状和变密度的三度体，本文提出了基于剖分算法的重力建模方法采用剖分算法将三维目标地质体分解为若干变密度四面体体元，推导了基于四面体体元的重力正演公式；以规则形状且连续密度变化的地质体（如长方体）和非规则形状变密度的地质体（如倾斜台阶组合模型体）为例，比较分析了常规块体方法和剖分算法应用于重力正演的计算效果，并采用共轭梯度法加密度约束条件对非规则形状变密度的地质体进行了反演计算模拟结果表明对于规则形状且连续密度变化的地质体，基于常规块体方法和剖分算法的重力建模均能获得满意的计算结果；对于非规则形状变密度的地质体，基于剖分算法的重力建模能获得更加合理的计算结果；随着对地质体分块数的增加，尽管常规块体方法的建模结果将逐渐逼近真实地质体